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SECTION  I 


T/tf  AtV  Foro*  .  _  INTRODUCTION 

The  U.S.  Air  Force  (USAF^Ls  taking  steps  to  improve  and  update 
^k6  safety  procedures  available  for  use  during  operations  involving 

hazardous  chemicals  and  fuels.  In  support  of  this  effort,  the  r  , 

i  jjo  -r 

Lawrence  Livermore  National  Labft»atory  < LLNL) >c onduc t e d^ja— eer ic s -of/ 

\tXtOiC  C.  {i~CLC!4c  3  le)  o£ 

pSltrogen  tetroxide  (N^O^)  e*  fur  thn  the^U.S. 

Department  of  Energy  (DOE)'-Nevada  Test  Site ^NTSHf?  during  1983.  "~The 
^0^  liquid  was  released  directly  on  the  ground  through  a  circular 
distribution  piping  network  with  six  outlets .  Although  the  test 
site  surface  is  known  for  water  impermeability,  absorbed 

readily,  causing  the  surface  to  buckle  and  heave  up  several  inches. 
At  the  completion  of  each  test,  nitrogen  gas  was  used  to  purge  the 
piping  system  of  N^O^,  thereby,  enhancing  ^0^  vaporization.  An 
array  of  sensors  was  placed  25  meters  downwind  of  the  gas  source  in 
order  to  determine  the  mass  evolution  rate  of  the  released  gas. 
Determination  of  the  mass  evolution  rate  was  to  be  accomplished  by 
measurement  of  the  NgO^  concentration,  temperature,  and  velocity  as 
the  vapor  cloud  passed  through  the  array.  The  chemical  interaction 
of  ^0^  v*por  with  ambient  humidity  and  oxygen  is  discussed  in 
Section  II.  In  Section  III,  a  range  of  source  evolution  rates  is 
obtained,  based  on  the  measured  mass  rates  at  25  meters  and  the 
conclusions  of  Section  II.  Another  array  of  sensors  designed  to 
measure  NO2  concentration,  temperature,  and  velocity  was  placed  785 
meters  downwind  of  the  source  to  measure  the  concentration  and 
vertical  extent  of  the  gas  cloud  as  it  moved  downwind  after  release. 


Concentration  measurements  at  785  meters  are  adjusted  for  reaction 
effects  based  on  the  measured  mass  evolution  rate.  Two  of  these 
releases.  Eagle  3  and  6,  provided  sufficient  data  for  comparison  and 
assessment  of  atmospheric  gas  dispersion  modeling  techniques 
(Reference  1) ;  test  conditions  for  Eagle  3  and  6  are  summarised  in 
Table  1. 

In  his  final  report,  McRae  (Reference  1)  showed  that  standard 
passive  dispersion  modeling  techniques  Inadequately  described  the 
dispersion  processes  present  during  Eagle  3  and  6.  (N^O^  vapor  is 

heavier  than  air  for  the  same  temperature  and  pressure.)  In  fact, 
it  is  well  established  that  atmospheric  dispersion  of  heavier-than- 
air  gases  (HTAGes)  may  differ  importantly  from  passive  dispersion 
(trace  contaminant  dispersion) .  The  following  phases  are  used  to 
describe  HTAG  dispersion  in  the  atmosphere  (References  3  and  4) : 

--negative  buoyancy-dominated  dispersion 

--stably  stratified  shear  flow 

--passive  dispersion  due  to  atmospheric  turbulence 
An  interactive  computer  model  (DEGADIS)  which  accounts  for  these 
dispersion  phases  has  been  developed  for  the  U.S.  Coast  Guard 
(Reference  3) .  Although  application  of  DEGADIS  has  been  primarily 
directed  to  the  prediction  of  concentrations  in  the  range  of 
hydrocarbon  lower  flammability  limits  (1  to  5  volume  percent) ,  an 
objective  of  this  research  is  to  begin  to  examine  the  applicability 
of  DEGADIS  to  the  prediction  of  concentrations  for  toxic  gases  by 
comparison  of  DEGADIS  predictions  with  the  Eagle  series  tests. 
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Eagle  3  Eagle  6 


Air  temperature  at  12  m  (°C) 

21.9 

22.6 

Wind  speed  at  12  m  (m/s) 

3.66 

5.58 

Wind  variability  at  12  m 

7.6° 

10.8° 

Barometric  pressure  (atm) 

0.9079 

0.9093 

Relative  humidity 

45% 

35% 

Reported  friction  velocity  (m/s) 

0.081 

0.148 

Total  mass  ^0^  released  (kg) 

6090 

4930 

Section  IV  compares  model  predictions  to  adjusted  concentrations 
and  to  observed  concentration  distributions  from  Eagle  3  and  6. 
Finally,  Section  V  examines  the  usefulness  of  a  release 
Richardson  number  for  determining  the  relative  importance 
of  the  three  dispersion  phases  mentioned  above. 
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SECTION  II 

CHEMICAL  INTERACTION  OF  N^/NOg  MIXTURES 
WITH  AMBIENT  HUMIDITY 


In  the  liquid  phase,  N20^/N02  mixtures  are  predominantly  N20^; 
at  1  atm,  the  concentration  of  N02  in  a  liquid  mixture  is  less  than 
0.2  mass  percent  (Reference  5).  In  the  vapor  phase,  the  equilibrium 
between  N20^  and  N02  expressed  as 


N2°4 


(8) 


2N0 


2  (g) 


(A) 


has  been  extensively  Investigated  (References  6 ,  7 ,  and  8) . 

Reaction  (A)  is  essentially  Instantaneous  (Reference  5) ,  and  the 

equilibrium  constant  K  .  based  on  ideal  gas  behavior  and  the 

eq,A 

thermodynamic  values  given  in  Table  2  Is  given  by 
2  2 

K  N°2  y”°2 
eq.A  P„  0  (1  atm)  “  y 

N2°4  N2°4 

-  exp  [23.48  -  0.3473  In  (T)  -  6983/T]  (1) 

where  p  is  the  partial  pressure  in  atm  and  y  is  the  mole  fraction;  T 
is  the  mixture  absolute  temperature  (K) . 
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TABLE  2.  THERMODYNAMIC  PROPERTIES  (Reference  9) 


AHf 

kcal/gmol 

ACf 

kcal/gmole 

C 

cal/gfiol  K 

NO 

(g) 

21.57 

20.69 

7.133 

NO. 

(g) 

7.93 

12.26 

8.89 

N2°4  ,  . 

<g) 

2.19 

23.38 

18.47 

HNO 

(g) 

-32.28 

-17.87 

12.75 

HNO 

(1) 

-41.40 

-19.10 

H2° 

(g) 

-57.796 

-54.634 

8.025 

A.  NEAR-FIELD  CHEMICAL  INTERACTIONS 

The  LLNL  infrared  (IR)  spectrophotometric  sensors  used  for  the 
25-meter  array  in  the  Eagle  series  tests  measured  the  transmission 
through  the  gas  sample  over  four  different  frequency  bands  denoted 
as  Channels  M,  F,  R,  and  E;  Channels  E  and  F  are  reference  channels 
because  no  absorption  at  these  frequencies  would  be  caused  by  the 
presence  of  any  of  the  gas  species  under  consideration.  Channels  M 
and  R  were  intended  to  measure  the  N^O^/NOg  mixture  concentration 
based  upon  the  additional  absorption  due  to  ^0^  (utilizing  the 
equilibrium  relationship  between  NOg  and  N^O^  given  above) . 

Unfortunately,  typical  IR  sensor  output  during  the  spill  tests 
indicated  essentially  the  same  attenuation  in  all  four  channels 
indicating  something  other  than  NgO^,  N02’  HN03’  HN®2'  an<*  amb*ent 
atmospheric  gases  was  present.  When  tested  using  NgO^  from  onsite 


storage,  the  IR  sensors  behaved  as  expected,  showing  little 
attenuation  in  Channels  E  and  F.  Grab  samples  of  the  gas  cloud  from 


the  25 -meter  array  and  in  the  spill  pipe  were  analyzed  at  LLNL 

by  mass  and  IR  spectroscopy;  no  indication  of  a  foreign  gas  capable 
of  producing  attenuation  in  all  four  channels  was  present.  As 
reported  in  Reference  2,  the  attenuation  was  concluded  to  be  the 
result  of  aerosol  scattering  which  would  produce  broad-band 
attenuation. 


A  possible  source  of  the  aerosol  is  from  the  reaction  of  ambient 
humidity  with  N20^/N02  to  form  nitric  acid  (HNO^)  which  has  been 
extensively  studied  (References  10,  11,  and  12).  Aerosol  formation 
is  known  to  take  place  when  the  partial  pressure  of  HNO^  exceeds  its 
vapor  pressure  at  the  mixture  temperature  (i.e.,the  mixture 


temperature  is  below  the  HNO^  dew  point  temperature) .  England  and 
Corcoran  (Reference  10)  proposed  the  following  mechanism  to  describe 
the  reaction  of  N204/N02  with  water  vapor  when  the  partial  pressure 
of  HNO^  does  not  exceed  its  vapor  pressure  (estimated  to  be  50  ppm 
of  HNOj  at  atmospheric  pressure  and  ambient  temperatures) : 


N2°4  (g)  -  2N02 


(A) 


N2°4  (g)  +  H2°  (g) 


«no2  (g)  +  hno3  (g) 


+  H2°  (g) 


Using  Table  2,  the  equilibrium  constant  without  mist  formation  is 
given  by 
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yHNO.  yN0 

Ke  D  “  ~3 — -  "  exp  t'12.95  -  1.038  In  T  +  4212/ i] 

q’  yN02  yH20  (2) 

When  the  concentration  of  HNO^  exceeds  50  ppm,  the  gas  phase 
reaction  supplies  nitric  acid  nuclei  where  water  vapor  readily 
condenses  to  form  strong  nitric  acid  droplets.  The  overall  reaction 
is  given  by 

3N0  +  H„0  -  NO  +  2HN0  (E) 

(8)  (g)  (g)  (aq> 

although  the  reaction  mechanism  remains  in  the  gas  phase 

(Reference  11).  As  Reactions  (D)  or  (E)  proceed,  the  ambient  water 

vapor  present  acts  to  reduce  the  N02  concentration,  thereby,, reducing 

the  quantities,  and  the  mist  formed  by  Reaction  (E)  further 

inhibits  the  measurement  of  N.O, . 

Z  4 

Although  the  mist  presence  made  direct  measurement  of  the 
total  cloud  concentration  impossible  with  the  sensors  which  were 
deployed,  McRae  et  al.  (Reference  2)  analyzed  the  IR  gas  sensor 
data  assuming  the  attenuation  was  the  same  in  channels  M  and  F . 
Under  this  assumption,  reported  concentrations  apply  to  the 
vapor  portion  of  the  cloud  only.  The  following  development 
attempts  to  provide  an  upper  and  lower  bound  of  the  total  cloud 
concentration  based  on  the  observed  concentration  of  the  vapor 
portion  of  the  cloud  assuming  no  reaction  takes  place.  These 
concentration  bounds  are  compared  with  the  vapor  dispersion 
models  in  Section  IV. 


fL->  <  1 If. 


;i 


If  the  gas  is  considered  to  be  a  mixture  of  air 
(nitrogen/oxygen  mixture)  xA,  water  (ambient  humidity)  xw,  and 
an  equilibrium  gas  phase  mixture  of  NO2/N2O4  xg  without  any 
reaction  taking  place,  the  sum  of  these  mass  fractions  would  be 

x  +  x  +  x  -  x  +  (1  +  w)  x.  -  1. 
s  w  A  s  A 


where  w  is  the  ambient  absolute  humidity.  Following  any  reaction, 
the  sum  of  the  mass  fractions  would  be 


x  +  x  +  x.,  +  x.  +  x  .  +  x  .  . 

w  A  w,A  acid 


where  the  mass  fraction  of  the  measured  mixture  is  x^,  the 

unabsorbed  and  unreacted  water  is  x^,  the  NO  present  due  to  reaction 

is  xN,  the  water  absorbed  at  the  HNO^  nuclei  is  xw  and  the  HNO^ 

present  due  to  reaction  is  x  ...  Equation  (4)  assumes  the  mass 

acid 

fraction  of  air  remains  unchanged  (open  to  question  for  longer  time 
scales  due  to  combination  of  0 2  with  NO,  and  discussed  later),  and 
there  is  little  error  in  the  measured  NOj/NjO^  concentration. 
Combining  Equations  (3)  and  (4), 


X  -  X  +x„  +  x  . ,  -  (x  -  X  -  x  .) 
s  v  N  acid  w  w  w,A 


where  (x  -  x  -  x  . )  Is  the  mass  of  water  reacted  per  mass  of 
w  w  w,  A 

mixture.  If  the  stoichiometry  of  Reactions  (D)  or  (E)  holds,  the 
following  ratios  apply: 

kg  NO  formed  _  1  kmol  NO  formed  30  kg  NO  formed 

XN  ”  kg  mixture  kmol  H^O  reacted  kmol  NO  formed 


i  kmol  H^O  reacted 
18  kg  HgO  reacted 


[1  fk®  H2° 

kg  m 

-  . 


reacted 


mixture 


* »*  •  >  V  ’  1  ' 

•  a  - mlT  ^  */> 


and 


kg  HNOj  formed 
xacid  "  kg  mixture 


Using  these  ratios. 


’2  kmol  HNO^  formed 

r63  kg  HN03  formed' 

kmol  H_0  reacted 

L  ^ 

kmol  HNO  formed 

1  kmol  H^O  reacted 

kg  H^O  reacted 

18  kg  H„0  reacted 

kg  mixture 

x 

s 


x  +  7.67  (x  -  x 
v  w  w 


Xv.A> 


(6) 


Although  Reaction  (E)  has  been  shown  to  be  fast  (References  11 
and  12),  the  reaction  has  not  been  shown  to  go  to  completion.  If  a 


conversion  efficiency  f  is  defined  as 

kmol  NC>2  reacted 
f  *  kmol  NO^  without  reaction 


where  the  (kmol  N02  without  reaction)  is  the  quantity  of  N02/N2°4 
without  any  reaction  taking  place  and  expressed  as  pure  NO^ ,  the 
mass  of  water  reacted  per  mass  of  mixture  is  given  by: 


(x  -  x  -  x  .) 

w  w  w,A 


kg  HgO  reacted 
kg  mixture 


18  kg  H^O  reacted 

kmol  H20  reacted 

kmol  H„0  reacted 

»  ^  - 

3  kmol  N0„  reacted 
^  * 

kmol  NOg  reacted 


kmol  NOj  without  reactionj 

kmol  NOj  without  reaction 
kg  mixture 


6f 


kmol  NC^  without  reaction 
kg  NOg/NjO^  without  reaction 

fkg  NOg/NgO^  without  reaction 


kg  mixture 


But,  the  second  term  in  brackets  is  x  while  the  first  term  is 

s 


(1/46).  Or, 


(7) 


(x  -  X  -  X 

w  w  w,A  46  s 


which  with  Equation  (6)  gives 


X  “  X 

s 


,/<!  "  f> 


Although  Equation  (7)  becomes  indeterminant  when  f  «  1,  this  does 
not  indicate  a  problem  since  if  f  *  1,  x  would  be  0.  It  should  be 
noted  that  f  will  be  less  than  1  if  water  is  the  limiting  reactant 
in  Reaction  (E) .  Furthermore,  no  assumption  about  which  reactant 
may  be  limiting  has  been  used  to  derive  Equation  (7). 

Goyer  (Reference  11)  measured  the  extent  of  Reaction  (E)  by 
mixing  humid  air  with  a  dry  nitrogen  stream  which  had  been  bubbled 
through  liquid  ^0^;  NOg  gas  concentrations  were  measured 
photometrically,  and  a  conversion  efficiency  r)  was  defined  as 

CN02i  '  CN02f 
1  ”  - n - 

"V 

3 

where  C„_  .  is  the  initial  NO.  concentration  (kmol  NO./m  )  and 
Nt>2i  t  i 

C„„  is  the  final  NO.  concentration.  Using  Goyer's  data  (Reference 
NO-f  l 


11),  t)  <  0.07  for  the  conditions  of  Eagle  3  and  Eagle  6.  It  can  be 
shown  that 

f  -  17/(1.  +  2yNo2l/Keq  A) 

where  y^  ^  is  the  initial  mole  fraction  of  N02>  Since  f  <  r\ 

<  0.07,  use  f  =  0.07  as  an  upper  bound  of  the  extent  of  reaction. 
Or, 
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depending  on  the  actual  value  of  f. 


B.  FAR- FIELD  CHEMICAL  INTERACTIONS 

In  addition  to  the  previously  mentioned  reactions,  NO  may  be 
oxidized  to  NOg  by  the  reaction 

2N0  +  02  •+  2N02  (F) 

Although  Reaction  (F)  is  suspected  to  be  unimportant  at  the  25 -meter 
sensor  array,  it  may  be  important  for  the  785-meter  sensor  array 
since  it  will  theoretically  go  to  completion.  Burdick  (Reference 
13)  investigated  Reaction  (F)  and  found  the  kinetics  to  be  described 


by 


dC 


NO 

dt 


Ki  \c*o 


with  the  value  of  ■  1,2  x  10^  1^/gmol^s  at  25°C.  Using  the 
maximum  possible  concentration  of  N02/N20^  at  the  25-meter  array  and 
assuming  the  concentration  of  oxygen  is  essentially  constant, 
the  minimum  value  of  (CM_/C1._  )  can  be  estimated  to  be  about  0.40 

NO 

for  Eagle  3  and  0.52  for  Eagle  6.  Although  Reaction  (F)  generates 
N02  at  the  expense  of  NO,  this  will  force  more  HNO^  to  be  generated 
due  to  the  equilibrium  constraints  of  the  faster  Reaction  (D)  (at 
the  expense  of  N02> .  The  net  effect  of  Reaction  (F)  is  thought  to 
be  small. 

Although  Goyer  (Reference  11)  observed  the  rapid  evaporation  of 
nitric  acid  mist  formed  for  relative  humidities  below  about  20 
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percent,  the  evaporation  process  is  believed  to  be  much  slower  due 
to  the  higher  humidity  of  Eagle  3  and  6.  It  should  be  noted  that 
the  first  HNO^  damage  to  the  low  level  IR  sensors  at  the  25-meter 
array  occurred  after  Eagle  3;  Eagle  3  was  the  first  test  with  a 
relative  humidity  much  greater  than  20  percent.  Since  Goyer's 
experiments  (Reference  11)  were  carried  out  with  air,  any 
decomposition  of  HNO^  to  NO^  would  be  taken  into  account  with  the 
observed  reaction  extent  (including  any  decomposition  acceleration 
due  to  the  presence  of  NO).  For  the  data  of  the  785-meter  array, 
the  same  mass  fraction  range  in  Equation  (8)  is  assumed  to  apply  to 
each  of  the  measured  concentrations. 


SECTION  III 


ADJUSTED  MASS  FLUX  AND  CONCENTRATION 

The  array  of  sensors  placed  25  meters  downwind  of  the  gas  source 
was  designed  to  measure  the  concentration,  temperature,  and 

velocity  of  the  gas  cloud  as  it  moved  downwind.  From  these 
measurements,  the  mass  rate  of  gas  (m)  25  meters  downwind  of  the 
source  is  given  by 

■  -  Ja  pXsUdA 

where  p  is  the  mixture  density  and  u  is  the  measured  velocity.  With 
a  steady-state  release,  the  mass  rate  of  gas  passing  a  plane  25 
meters  downwind  of  the  source  is  also  the  mass  evolution  rate  of  gas 
from  the  source.  Note  that  the  mass  flux  of  gas  is  based  on 
NO^/^O^  concentrations  without  any  reaction  taking  place. 

The  concentration  range  from  Equation  (8)  can  be  used  to  adjust 
the  measured  vapor  flux  at  the  25-meter  array.  Using  the  (Reference 
1)  measured  vapor  evolution  rate  of  NO^/NgO^  for  Eagle  3  and  6 
(McRae,  1985)  and  Equation  (8),  the  mass  evolution  rate  for  Eagle  3 
ranged  between  2.9  and  3.1  kg/s  (170  kg/min  and  190  kg/min)  while 
the  mass  evolution  rate  for  Eagle  6  ranged  between  1.6  and  1.7  kg/s 
(97  kg/min  and  103  kg/min) .  These  ranges  of  the  source  evolution 
rate  are  consistent  with  the  observed  extensive  outgassing  times  for 
both  Eagle  3  and  6.  McRae  (Reference  1)  reported  that  100  to  1000 
gallons  of  water  were  sprayed  on  the  release  area  to  dilute  the  NgO^ 
absorbed  into  the  ground  so  that  acceptable  vapor  concentration 
levels  would  be  realized.  Based  on  these  evolution  rates,  78  to  80 
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percent  of  the  released  for  both  Eagle  3  and  Eagle  6  absorbed 

into  the  soil  at  the  source  during  the  test,  in  agreement  with 
previous  small  scale  tests  at  Edwards  Air  Force  Base  (Reference  14) . 
(This  agreement  may  not  be  significant  due  to  the  differences  in 
the  soil  between  the  Eagle  and  Edwards  tests.) 

In  addition  to  the  measured  MO^  concentration  at  785  meters 
shown  in  Table  3,  the  following  method  was  used  as  another  estimate 
of  the  N(>2  concentration  at  785  meters.  For  each  of  the  Eagle  tests 
reported  by  McRae  (Reference  1) ,  an  appropriate  steady  source  time 
period  was  chosen  for  each  experiment,  taking  into  account  ambient 
wind  velocities  and  directions  over  the  travel  time  between  the 
source  and  the  785 -meter  array.  After  choosing  this  time,  the 
observed  concentration  contours  were  adjusted  to  a  Gaussian  equiva¬ 
lent  cross  section.  These  cross  sections  were  integrated  to  deter¬ 
mine  the  observed  mass  flux  of  NO^/NjO^  through  the  array.  Since  a 
steady- state  experiment  must  observe  the  same  mass  flux  between  the 
25-meter  array  and  the  785-meter  array,  a  maximum  centerline 
concentration  may  be  calculated  assuming  the  observed  profiles  (a^ 
and  0^  values)  to  be  applicable.  The  results  of  these  calculations 
are  summarized  in  Table  3.  All  of  the  reported  concentrations  are 

in  equivalent  concentration  of  N0„  (yNn  +  2y  ) ,  and  for  these 

l  NU2  n2°4 

concentration  ranges  at  785  meters,  y^Q  »  yN  ^  The  overall 
N02  concentration  range  in  Table  3  represents  the  minimum  and 
maximum  N02  concentrations  from  all  indications  for  the  steady  state 
N02  concentration  at  785  meters  for  Eagle  3  and  6. 
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TABLE  3.  SUMMARY  OF  SOURCE  RATES  AND  DOWNWIND  CONCENTRATIONS  FOR 
EAGLE  3  AND  6 


Eagle  3  Eagle  6 


Mass  rate  range  of  unreacted 

N02>,li204  at  ^  meters  (kg/s)  2. 9-3.1  1.6-1. 7 

Total  mass  of  NO./N.O,  passing 

25-meter  array  (kg)  1200-1300  970-1030 

Observed  Gaussian  equivalent 
distribution  at  785  meters 
a  (crosswind)* 

1 7 ^  (vertical)* 

Observed  maximum  NO.  concentration 

at  785  meters  (ppm)*  500  315 

Observed  maximum  NO^  concentration  at 
785  meters  taking  reaction  into 

account  (ppm)  535  340 

Range  of  maximum  NO^  concentration  at 
785  meters  based  on  mass  flux  at 

25  meters  (ppm)  970-1040  160-170 

Overall  NO^  concentration  range 
at  785  meters  (ppm) 


35  m  35  m 

3.8m  7.6m 


500-1040 


160-340 


SECTION  IV 


GAS  DISPERSION  MODEL  COMPARISONS 

McRae  (Reference  1)  compared  several  Gausslan-prof lie  model 
predictions  to  the  measured  downwind  NO^  concentration  as  well 
as  the  equivalent  Gaussian  cross-section  profiles  observed  in 
the  experiments.  Two  passive  dispersion  models  (Ocean 
Breeze/Dry  Gulch  and  Pasquill -Hanna)  are  compared  to  the  data 
presented  in  Section  III.  Also,  comparisons  are  made  to  the 
DEGADIS  model  developed  for  the  U.S.  Coast  Guard  (Reference  3). 

The  Ocean  Breeze/Dry  Gulch  (OB/DG)  model  is  the  product  of 
two  experimental  programs  involving  the  release  and  measurement 
of  a  trace  atmospheric  contaminant.  The  results  were  correlated 
to  a  simple  empirical  model  (Reference  15) .  The  same  procedures 
used  by  McRae  (Reference  1)  were  used  for  calculating  the  wind 
field  variation  {a ^)  and  the  characteristic  temperature 
difference  (AT) .  Although  the  OB/DG  model  predicts  the  - 
concentration  at  1.5  meters  elevation,  the  difference  between 
this  elevation  and  the  sensor  elevation  of  1.0  meters  would  not 
be  great.  The  OB/DG  results  for  Eagle  3  and  6  are  listed  in 


Table  4. 
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TABLE  4.  COMPARISON  OF  EAGLE  3  AND  EAGLE  6  TEST  RESULTS  AND  GAS 
DISPERSION  MODEL  PREDICTIONS 


Maximum 

NO.  Concentration  ay 

Range  (ppm)*  (m) 

oz 

(m) 

EAGLE  3 

Test  Results 

500-1040 

35 

3.8 

|  OB/DG 

68-73 

-- 

-- 

"  Gaussian  Plume 

1* 

<• 

68-73 

60.5 

31.9 

DEGADIS 

M 

880-1170 

57.6-60.7** 

2. 3-2. 9*** 

,1 

jj  EAGLE  6 

S  Test  Results 

160-340 

35 

7.6 

> 

£  OB/DG 

20-21 

-- 

-- 

1  Gaussian  Plume 

25-27 

60.5 

31.9 

5  DEGADIS 

ft 

% 

190-220 

55.6-56.6** 

4. 5-4. 9*** 

*The  concentration  range  for  the  model  predictions  are  for  the 
estimated  source  evolution  rate  range. 


**calculated  as  S 


★♦★calculated  as  S  /Jl 
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The  Pasquill -Hanna  Gaussian  plume  model  has  proven  to  be 
applicable  to  atmospheric  dispersion  problems  when  the 
dispersion  of  the  contaminant  is  only  a  function  of  the 
atmospheric  turbulence  and  the  plume  does  not  perturb  the 
ambient  flow  field.  The  steady  Gaussian  plume  model  for  ground 
level  releases  is  given  by 


iVlV  . .  l  •  .V.'-.V 
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ss 
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y  - - exp 

v  jrp  a  a  u  r 
v  y  z 


where  y  is  the  contaminant  mole  fraction,  m  is  the  source 

•'v 

evolution  rate,  p is  the  vapor  density  at  ambient  temperature 
and  pressure,  u  is  the  average  wind  speed,  and  a and 
determine  the  shape  of  the  Gaussian  profile  taken  from  Hanna  et 
al.  (Reference  16).  The  same  procedures  used  by  McRae 
(Reference  1)  were  used  in  determining  the  parameters  in 
Equation  (10) .  The  model  results  for  Eagle  3  and  6  are  listed 
in  Table  4. 

The  DEGADIS  (DEnse  GAs  Dispersion)  model  is  an  adaptation  of 
the  Shell  HEGADAS  model  described  by  Colenbrander  (Reference  17) 
and  Colenbrander  and  Puttock  (Reference  18).  The  near-field, 
buoyancy- dominated  regime  is  modeled  using  a  lumped  parameter 
model  of  a  heavier-than-air  gas  (HTAG)  "secondary  source"  cloud 
which  incorporates  air  entrainment  at  the  gravity -spreading 
front  using  a  frontal  entrainment  velocity.  The  downwind 
dispersion  phase  of  the  calculation  assumes  a  power  law 
concentration  distribution  in  the  vertical  direction  and  a 
modified  Gaussian  profile  in  the  horizontal  direction  with  a 
power  law  specification  for  the  wind  profile  (Figure  1).  A 
description  of  the  DEGADIS  model  is  included  in  Appendix  A.  The 
source  model  represents  a  spatially  averaged  concentration  of 
gas  present  over  the  primary  source,  while  the  downwind 
dispersion  phase  of  the  calculation  models  an  ensemble  average 


of  the  concentration  downwind  of  the  source.  The  vertical 
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Figure  1.  Schematic  Diagram  of  DEGADIS  Dense  Gas 
Dispersion  Model 


mixing  rate  is  based  on  laboratory* scale  data  for  vertical 
mixing  in  stably  density- stratified  fluids  reported  by  Kantha  et 
al.  (Reference  19),  Lofquist  (Reference  20),  and  McQuaid 


(Reference  21) .  The  vertical  dispersion  parameter  S  and  the 

z 

horizontal  dispersion  parameter  determine  the  vertical  and 
horizontal  profiles,  respectively.  The  rates  of  change  of 
and  approach  the  rates  of  change  of  a and  a respectively, 
as  the  density  of  the  plume  approaches  the  ambient  air  density. 

Because  several  parameters  are  calculated  from  the  Pasquill 
stability  class  and  the  surface  roughness,  DEGADIS  requires  no 
more  atmospheric  information  than  the  Pasquill-Hanna  Gaussian 
plume  model.  However,  because  of  the  dependence  of  the  vertical 
dispersion  rate  on  the  gas  mixture  density,  the  gas  mixture 
density  must  be  specified  as  a  function  of  concentration.  For 
locations  near  the  source,  the  gas  mixture  is  predominantly 
N20^,  while  for  distances  far  from  the  source,  the  gas  mixture 
is  predominantly  NO^ .  Two  simulations  were  made  to  bound  the 
behavior  of  the  reacting  gas  mixture  as  it  moved  downwind.  For 
one  bounding  simulation,  the  contaminant  gas  was  assumed  to  be 
pure  NOg  at  the  temperature  of  the  liquid  pool 
(~7°C  for  Eagle  3  and  15°C  for  Eagle  6).  For  the  other  bounding 
simulation,  the  contaminant  gas  was  assumed  to  be  pure  at 

the  temperature  of  the  liquid  pool;  under  these  conditions,  the 
predicted  mole  fraction  was  multiplied  by  2  for  the 

equivalent  N02  mole  fraction.  For  these  bounding  simulations, 


the  calculations  were  performed  assuming  no  reaction  took  place 


,1 

I 


(tinea  any  reaction  which  would  occur  would  cause  the  gas 

mixture  density  to  fall  between  the  bounds  of  without 

reaction  and  without  reaction).  The  DEGADIS  simulation 

results  for  Eagle  3  and  6  are  presented  in  Table  4.  The  ranges 

of  NO.  concentration,  a  ,  and  a  represent  the  minimum  and 
€•  y  z 

maximum  variable  values  calculated  using  the  source  evolution 

rate  range  and  both  of  the  bounding  density  specifications. 

For  Eagle  3  and  6,  the  passive  dispersion  models  tend  to 

underpredict  the  maximum  concentration  because  of  an 

overprediction  of  the  vertical  mixing  rate  (i.e.,the  predicted 

a  is  greater  than  the  observed  a  ) .  On  the  other  hand,  DEGADIS 
z  z 

predicts  maximum  centerline  concentrations  and  similarity 
distribution  parameters  a  and  a  which  are  consistent  with  the 

y  * 

experimental  data;  closer  agreement  of  DEGADIS-predicted  values 
of  az  with  observed  values  of  a g  compared  to  the  passive 
Gaussian  plume  models  is  attributed  to  correctly  predicting  the 
decreased  vertical  mixing  rate  for  the  initial  phase  of  the 


release. 


SECTION  V 


RELEASE  RICHARDSON  NUMBER 

The  fact  that  passive  dispersion  prediction  techniques  may 
be  inadequate  for  releases  of  N2O4  might  be  anticipated  because 
N2O4  vapor  is  heavier  than  air  for  the  same  temperature  and 
pressure.  Passive  dispersion  prediction  techniques  are  based  on 
the  assumption  that  the  dispersing  contaminant  does  not  disturb 
the  ambient  flow  field.  In  contrast,  the  dispersion  phases  of  a 
heavier-than-air  gas  (HTAG)  may  also  include  negative  buoyancy- 
dominated  dispersion  and  stably  stratified  shear  flow  as  well  as 
passive  dispersion  (References  3  and  4).  The  negative  buoyancy- 
dominated  dispersion  phase  is  characterized  by  turbulent  flow 
structures  (for  example,  the  classic  "gravity  current")  which 
result  in  rapid  mixing.  Between  these  two  extremes,  the  stably 
stratified  shear  flow  phase  is  characterized  by  vertical  mixing 
rates  which  are  lower  than  the  passive  dispersion  phase,  and  the 
turbulent  structures  in  the  negative  buoyancy- dominated 
dispersion  phase  are  generally  no  longer  significant. 

Notice,  however,  that  these  dispersion  phases  are  not 
governed  only  by  the  density  of  the  dispersing  gas  but  also  by 
other  conditions.  Consider  the  release  of  a  gas  slightly 
heavier  than  air  in  a  very  light  wind;  such  a  release  could 
exhibit  the  characteristics  of  the  negative  buoyancy- dominated 
dispersion  phase  even  though  the  contaminant  is  only  slightly 
heavier  than  air.  Furthermore,  consider  a  small  release  of  gas 
which  is  several  times  heavier  than  air  in  high  winds;  such  a 
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release  could  be  controlled  by  the  ambient  wind  field  (passive 
dispersion)  even  though  the  contaminant  is  several  times  heavier 
than  air.  Therefore,  the  question  of  whether  a  release  is 
dominated  by  "dense  gas  effects"  is  governed  by  three  factors, 
including  the  density  of  the  released  gas,  the  rate  of  release, 
and  the  ambient  flow  field.  (A  release  is  said  to  be  dominated 
by  "dense  gas  effects"  if  the  maximum  concentration  or  the 
concentration  distribution  at  the  concentration  of  interest  is 
significantly  affected  by  the  negative  buoyancy-dominated 
dispersion  phase  and/or  the  stably  stratified  shear  flow  phase.) 

These  qualitative  characteristics  can  be  quantitatively 
predicted  by  the  release  Richardson  number 


where  p is  the  initial  (source)  gas  density.  This  Richardson 
number  represents  a  ratio  of  the  characteristic  potential  energy 
of  the  gas  release  to  a  measure  of  the  ambient  turbulent  kinetic 
energy.  Vhen  the  characteristic  potential  energy  of  the  release 
is  much  greater  than  the  ambient  turbulent  kinetic  energy 
(Ric  »  1),  the  negative  buoyancy- dominated  dispersion  phase 
predominates  near  the  source.  As  the  gas  moves  downwind  and 
dilutes,  the  dispersion  is  next  controlled  by  the  stably 
stratified  shear  flow  phase,  and  as  the  gas  is  further  diluted, 
the  dispersion  is  then  passive.  In  contrast,  when  the 
characteristic  potential  energy  of  the  gas  release  is  much  less 


than  the  ambient  turbulent  kinetic  energy  (Ric  <  1) ,  the  passive 
dispersion  phase  predominates  near  the  source.  Based  on 
experimental  work  reported  by  Britter  (References  22  and  23) , 
three  release  Richardson  number  ranges  roughly  correspond  to  the 
three  dispersion  phases  discussed  above,  namely: 

Ric  >  32  negative  buoyancy-dominated  dispersion 

1  <  Ric  <  32  stably  stratified  shear  flow 

Ri^  <  1  passive  dispersion 

where  the  intermediate  range  is  a  transition  zone  between  the 
two  extremes . 

For  Eagle  3  and  6,  p was  calculated, based  on  an  equilibrium 

mixture  of  NO^/N^O^  at  the  liquid  pool  temperature.  Using  this 

definition,  Ric  -  40  for  Eagle  3  and  Ric  -  4  for  Eagle  6.  (The 

Eagle  6  Ri  is  lower  than  the  Eagle  3  Ri  because  Eagle  6  had  a 
c  c 

lower  release  rate  and  a  higher  wind  speed  than  Eagle  6 . )  From 
these  Richardson  numbers,  the  dispersion  in  Eagle  3  and  6  would 
be  expected  to  be  Inadequately  described  by  passive  dispersion 
prediction  techniques  such  as  OB/DG  and  the  Pasquill -Hanna 
Gaussian  plume  model;  this  conclusion  is  consistent  with 
comparisons  of  the  range  of  observed  maximum  concentrations  and 
the  similarity  distribution  parameters  a  and  a ^  with  the  OB/DG 
and  Pasquill -Hanna  Gaussian  plume  model  predictions.  In 
contrast,  since  DEGADIS  was  designed  to  include  the  three 
dispersion  phases  of  a  HTAG,  the  dispersion  present  in  Eagle  3 
and  6  would  be  expected  to  be  adequately  described  by  DEGADIS. 
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Indeed,  as  indicated  in  Section  IV,  the  DEGADIS  predictions  are 
consistent  with  the  range  of  observed  maximum  concentrations  and 
the  parameters  Oy  and  az  for  Eagle  3  and  6.  Reference  3 
contains  further  comparison  of  DEGADIS  vith  other  field  tests. 
(If  Ric  <  1,  the  predictions  of  the  Pasquill-Hanna  Gaussian  plume 
model  and  DEGADIS  would  be  expected  to  be  consistent  between 
themselves  as  well  as  with  field-scale  data  (Reference  3).) 


SECTION  VI 


CONCLUSIONS 

Two  field-scale  releases  of  NjO^  (Eagle  3  and  6)  performed 
by  the  Lawrence  Livermore  National  Laboratory  for  the  U.S.  Air 
Force  as  part  of  the  Eagle  series  were  examined. 

An  analysis  of  the  chemical  interaction  of  N2O4  with  ambient 
humidity  and  oxygen  was  made.  Based  on  this  analysis,  the 
source  mass  evolution  rate  m  for  Eagle  3  was  estimated  to  be  2.9 
kg/s  <  m  <  3.1  kg/s  and  1 . 6  kg/s  <  m  <  1.7  kg/s  for  Eagle  6 . 

Based  on  these  mass  evolution  rates,  approximately  80  percent  of  the 
N2O4  released  in  Eagle  3  and  6  was  absorbed  into  the  ground  at 
the  spill  site. 

Observed  downwind  concentrations  and  Gaussian  equivalent 

concentration  profiles  ( a ^  and  a were  compared  to  model 

predictions  using  the  Ocean  Breeze/Dry  Gulch  model  (OB/DG) ,  the 

Fasquill -Hanna  Gaussian  plume  model,  and  DEGADIS.  OB/DG  and  the 

Pasquill -Hanna  model  generally  underpredicted  the  concentration 

due  to  an  overprediction  of  the  vertical  mixing  present.  On  the 

other  hand,  DEGADIS  predictions  were  consistent  with  observed 

concentrations  and  profile  parameters  a  and  a  ,  indicating  the 

y  z 

importance  of  density  stratification  in  these  tests. 

The  importance  of  density  stratification  In  atmospheric 

dispersion  can  be  predicted  by  analyzing  a  release  Richardson 

number  RI  .  A  heavier-than-air  gas  (HTAG)  release  is  dominated 
c 

by  "dense  gas  effects"  for  Ri^  >  1,  but  the  release  will  be 
passive  from  the  source  for  Ric  <1.  (A  HTAG  release  is  said  to 
be  dominated  by  dense  gas  effects  if  the  maximum  concentration 
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or  the  concentration  distribution  at  the  concentration  of 


interest  is  significantly  affected  by  the  negative  buoyancy- 
dominated  dispersion  phase  or  the  stably  stratified  shear  flow 
phase.)  For  Eagle  3  and  6,  Ri^  >  1  which  indicates  that  passive 
dispersion  prediction  techniques  could  Inadequately  predict  the 
dispersion  processes  present  in  these  tests;  this  conclusion  was 
substantiated  by  comparison  of  test  results  with  0B/D6  and  the 
Pasquill -Hanna  model  as  mentioned  above.  In  contrast,  DEGADIS, 
which  accounts  for  the  negative  buoyancy-dominated  phase  and  the 
stably  stratified  shear  flow  phase  as  well  as  the  passive 
dispersion  phase,  predicted  maximum  concentrations  and 
distribution  parameters  which  were  consistent  with  observed 
values . 
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APPENDIX  A 


DESCRIPTION  OF  THE 
DEGADIS  DENSE  GAS  DISPERSION  MODEL 


This  appendix  is  a  self-contained  document  with  its  own 
format  as  well  as  internally  consistent  numbering  system 
for  equations  and  figures. 


>  *—  U'fiVJjiJiy. 


A' >»■  -*n  .V  .fa  .*i  ,*»  -*■  ■>«  .s  >'..*<■  .*■  .»**  .**  ■» 


DESCRIPTION  OF  THE 
DEGADIS  DENSE  GAS  DISPERSION  MODEL 

The  DEGADIS  (DEnse  GAs  Dispersion)  model  was  developed  from 
research  sponsored  by  the  U.S.  Coast  Guard  and  the  Gas  Research 
Institute  (Reference  3).  DEGADIS  is  an  adaptation  of  the  Shell 
HEGADAS  model  described  by  Colenbrander  (Reference  17)  and 
Colenbrander  and  Puttock  (Reference  18).  DEGADIS  also  incorporates 
some  techniques  used  by  van  Ulden  (Reference  23) . 

If  the  primary  source  (gas)  release  rate  exceeds  the  maximum 
atmospheric  takeup  rate,  a  heavier -than- air  (HTAG)  blanket  is  formed 
over  the  primary  source.  This  near-field,  buoyancy- dominated  regime 
is  modeled  using  a  lumped  parameter  model  of  a  HTAG  "secondary 
source"  cloud  which  incorporates  air  entrainment  at  the  gravity- 
spreading  front  using  a  frontal  entrainment  velocity.  If  the 
primary  source  release  rate  does  not  exceed  the  maximum  atmospheric 
takeup  rate,  the  released  gas  is  taken  up  directly  by  the  atmosphere 
and  dispersed  downwind.  For  either  source  condition,  the  downwind 
dispersion  phase  of  the  calculation  assumes  a  power  law  concen¬ 
tration  distribution  in  the  vertical  direction  and  a  modified 
Gaussian  profile  in  the  horizontal  direction  with  a  power  law 
specification  for  the  wind  profile  (Figure  A-l).  The  source  model 
represents  a  spatially  averaged  concentration  of  gas  present  over 
the  primary  source,  while  the  downwind  dispersion  phase  of  the 
calculation  models  an  ensemble  average  of  the  concentration  downwind 
of  the  source . 

A.  HEAVIER- THAN -AIR  GAS  (HTAG)  SOURCE  CLOUD  FORMATION 

A  lumped  parameter  model  of  the  formation  of  the  HTAG  source 
cloud  or  blanket,  which  may  be  formed  from  a  primary  source  such  as 
an  evaporating  liquid  pool  or  otherwise  specified  ground  level 
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Figure  A-l,  Schematic  diagram  of  DEGADIS  Dense  Gas 
Dispersion  Model 


emission  source,  or  by  an  initially  specified  gas  volume  of 
prescribed  dimensions  for  an  instantaneous  release,  is  illustrated 
in  Figure  A-l.  The  gas  blanket  is  represented  as  a  cylindrical  gas 
volume  which  spreads  laterally  as  a  density-driven  flow  with 
entrainment  from  the  top  of  the  source  blanket  by  wind  shear  and  air 
entrainment  into  the  advancing  front  edge.  The  source  blanket  will 
continue  to  grow  over  the  primary  source  until  the  atmospheric 
takeup  rate  from  the  top  is  matched  by  the  air  entrainment  rate  from 
the  side  and,  if  applicable,  by  the  rate  of  gas  addition  from  under 
the  blanket.  Of  course,  the  blanket  is  not  formed  if  the 
atmospheric  takeup  rate  is  greater  than  the  evolution  rate  of  the 
primary  source.  For  application  of  the  downwind  calculation 
procedure,  the  blanket  is  modeled  as  being  stationary  over  the 
center  of  the  source  (x  -  0) . 

1.  SECONDARY  SOURCE  BLANKET  EXTENT  FOR  GROUND  LEVEL  RELEASES 

If  a  HTAG  blanket  is  present,  the  (downwind)  emission  rate  from 

the  blanket  is  equal  to  the  maximum  atmospheric  takeup  rate.  That 
2 

is,  for  E(t)/jrR  (t)  >  Q  ,  a  source  blanket  is  formed  over  the 
p  *max 

primary  source.  The  blanket  frontal  (spreading)  velocity  is  modeled 
as 

(A-l) 

where  p  is  the  average  density  of  the  source  blanket.  This  gravity 
intrusion  relationship  is  applicable  only  for  P  >  P&'>  the  value  of 
Cg  used  is  1.15  based  on  laboratory  measurements  of  cloud  spreading 
velocity  (Reference  3) . 

The  blanket  radius4R,as  a  function  of  time  is  determined  by 
integrating  dR/dt  -  u^.  For  ground  level  sources,  the  blanket 
spreading  is  stopped  (dR/dt  -  0)  when  the  total  mass  of  the  cloud  is 
decreasing  with  time.  The  radius  of  the  blanket  is  constrained  to 


be  greater  than  or  equal  to  the  radius  R^  of  any  primary  (liquid) 
source  present. 

2.  SECONDARY  SOURCE  BLANKET  EXTENT  FOR  INSTANTANEOUS  RELEASES 


The  gravity  intrusion  relationship  (Equation  (A-l))  will 
overpredict  initial  velocities  for  instantaneous,  aboveground 
releases  of  a  HTAG  since  no  initial  acceleration  phase  is  Included. 
In  this  case,  the  following  procedure  adapted  from  van  Ulden 
(Reference  23)  is  recommended. 

For  instantaneous  gas  releases ,  the  radially  symmetric  cloud  is 
considered  to  be  composed  of  a  tall  section  with  height  Hfc  and 
radius  R^  and  a  head  section  with  height  (Figure  A-2) .  A 
momentum 
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Figure  A-2.  Schematic  Diagram  of  a 

Radially  spreading  Cloud. 


balance  is  used  to  account  for  the  acceleration  of  the  cloud  from 
rest;  the  effect  of  ambient  (wind)  momentum  is  ignored.  Although 
the  following  equations  are  derived  assuming  the  primary  source  rate 
is  zero,  the  resulting  equations  are  assumed  to  model  the  secondary 
source  cloud  development  when  the  primary  source  rate  is  nonzero. 
When  the  frontal  velocity  from  the  momentum  balance  is  the  same  as 
Equation  (A-l),  the  momentum  balance  is  no  longer  applied  and  the 
frontal  velocity  is  given  by  Equation  (A-l). 

There  are  three  main  forces  acting  on  the  cloud:  a  static 
pressure  force  (F^) ,  a  dynamic  drag  force  (F^) ,  and  a  force  which 
accounts  for  the  acceleration  reaction  of  the  ambient  fluid, 
represented  as  a  rate  of  virtual  momentum  change  with  respect  to 
time  (-dP ^/dt) .  Denoting  the  momentum  of  the  head  and  tail  as 
and  P  respectively,  the  momentum  balance  is 


dP  =  d_ 
dt  dt 


(ph  ♦  Pt> 


F  +  F . 
p  d 


(A-2) 


or 


d_ 

dt 


+  P^  +  P  ) 
t  v 


F  +  F . 
p  d 


(A-3) 


The  terms  in  the  momentum  balance  are  evaluated  differently  for 
early  times  before  a  gravity  current  head  has  developed  (H^  <  Hfc) 
and  for  times  after  the  head  has  developed  but  the  cloud  is  still 
accelerating  (Figure  A-2).  Because  the  gravity  current  head 
develops  so  rapidly,  the  model  equations  describing  the  times  after 
the  gravity  current  head  forms  (H^  >  Ht)  are  derived  first.  The 
model  equations  describing  earlier  times  (H^  <  Ht)  use 
simplifications  of  the  equations  for  >  H  . 

a.  Unsteady  Gravity  Current 

When  the  cloud  accelerates  to  the  point  that  >  Ht  (Figures 
A-2,  A-3),  the  frontal  velocity  is  determined  from  the  momentum 
balance  (Equation  (A-2))  as  follows. 
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pressure  over  the  boundary  of  the  current.  Is 


-  [  \  gApHt  j  [  2*RHt  ]  -  *gApRH* 


(A-4) 


Neglecting  the  shear  stress  at  the  bottom,  the  dynamic  force  on  the 
current  Is  the  stun  of  the  drag  force  on  the  head  of  the  current  and 
the  lift  force  that  arises  due  to  asymmetry  in  the  ambient  flow 
around  the  head.  The  drag  force  is  represented  by 


FD  -  -  r  PA  ( 2*VA  )  "  -avVRVauf  (A*5 

where  d  is  an  effective  drag  coefficient  and  the  constant  a  is  an 
v  •  v 

empirical  ratio  of  the  average  head  depth  to  Hb  (av  -  W  • 

The  horizontal  acceleration  reaction  (-dP ^dt)  is  approximated 
by  the  reaction  to  an  accelerating  elliptical  cylinder  with  an 
aspect  ratio  H/R  (Reference  24) : 


(A-6) 


-  -  M  kio.’RK2"f ) 


and  the  vertical  acceleration  reaction  is  represented  as 


dp  ; 

v  d 

dt  dt 

z 

k2pa*RH2uf  ] 

(A-7) 

where  k^  and  k^  are 
constant.  Equations 

coefficients  of 

(A-6)  and  (A-7) 

order  one. 

give 

Using  a  single 

dP 

V 

-  tz —  **  -  e  itp 
dt  v  a 

d(RH2uf) 

(A- 8) 

dt 

Using  Equations  (A-4),  (A-5),  and  (A-8),  the  momentum  balance 
(Equation  (A-2))  becomes 

dP  2  2  d(RH2u  ) 

di  "  -  SVWf  *  VPa  - dt -  (A'9) 

Following  van  Ulden  (References  23,  25),  it  is  assumed  that  the 
potential  energy  decrease  due  to  slumping  of  the  cloud  is  offset  by 
the  production  of  kinetic  energy,  which  through  the  action  of  shear, 
is  partly  transformed  to  turbulent  kinetic  energy.  Part  of  the 
turbulent  kinetic  energy  is  transformed  back  into  potential  energy 
due  to  entrainment  of  air  by  the  cloud.  This  "buoyant  destruction" 
of  kinetic  energy  is  assumed  to  be  proportional  to  the  rate  of 
production  of  turbulent  kinetic  energy,  and  following  Simpson  and 
Britter  (Reference  26)  it  is  assumed  that  the  turbulent  kinetic 

3 

energy  production  rate  scales  as  vp  HRu_.  Then, 

&  £ 

|  gApH  -  e*paHRuJ  (A-10) 

which  can  be  written 
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^  .'-A 


dy  €  (2jrRH)uf  e(2*RH)uf 


(A-ll) 


where  e  is  an  empirically  determined  coefficient, 
represents  the  air  entrainment  rate, 


Noting  dV/dt 


—  «  € (2wRH)u. 
P  t 

j 


where  represents  the  air  entrainment  mass  rate. 


(A-12) 


i.J 


The  volume  integral 


I 


V  -  2jt  h(r ,t)rdr  (A-13) 

J0 

where  h(r,t)  is  to  be  expressed  in  terms  of  and  Ht>  and  the 
momentum  integral 


P  -  2?r  I  pu(r,t)h(r ,t)rdr  -  P_  +  P. 

J0  t  h 


(A-14) 


are  then  approximated  with  separate  analyses  of  the  head  and  tail  of 
the  current. 

In  the  tail  of  the  current,  the  shallow  water  equations  are 
assumed  applicable.  It  is  assumed  that  the  shape  of  the  current  is 
quasi- stationary  in  time,  and  the  layer -averaged  density  difference 
is  assumed  horizontally  uniform.  It  follows  that  the  volume  and 
momentum  of  the  tall  are  given  by 


vt  -  **J  ( +  “h  ] ' 2- 


(A-15) 


pt  '  5  -  (  5  "t  +  "h  ]  ‘I 


(A-16) 


uZ*: 


A  momentum  balance  for  the  head  region.  Figure  A-4,  assuming 
quasi-steady  state,  indicates  that  the  static  and  dynamic  pressure 
forces  on  the  head  should  be  balanced  by  the  net  flux  of  momentum 
due  to  flow  into  and  out  of  the  head.  The  static  pressure  and  drag 
are,  respectively 


fp  '  ( I  )  ( 2’Vh  )  - 
FD  -  -  dv  (  I  Vf  )  [  2"V  »vV  . 


(A-17 


"  '  WaWHi  <A-18 

Near  the  surface,  the  inward  flow  (u^  in  Figure  A-4)  carries 
momentum  into  the  head,  while  the  return  flow  (u^  in  Figure  A-4) 
carries  momentum  out  of  the  head.  Assuming  u^  “  u^,  =*  1/2  H^, 

and  u^  “  ^vuf*  t^ie  «“°™®ntum  flux  into  the  head  is  approximately 


A  [  2’8h"h 


(A-19 


H,  -f - /-  - 


-V-Hh-^U3 


Figure  A-4.  The  head  of  a  steady  gravity  current 
(References  26,  23). 


Upon  rearranging,  the  momentum  balance  on  the  head  gives 


-  l./f  da  -2  S2  ) 
gApH^  l  v  v  v  J 


-  C*  (A-20) 

when  5^  »■  0.2  and  d^  -  0.64;  Equation  (A-20)  then  specifies  the  head 
velocity  boundary  condition.  The  volume  of  the  head  is  determined 
by  assuming  that  the  head  length  scales  with  H^.  It  follows  that 

R  -  -  bv«1  (A- 21) 


where  b^  is  an  empirical  constant,  and  the  volume  of  the  head 
becomes 

\  -  '*vbv  <»  ♦  v  4  <*-22> 


If  the  layer -averaged  velocity  is  assumed  to  increase  linearly  with 
r,  it  follows  that 


%  "  Uf 


(A-23) 


and 


2n 

r 


-A 


R 


[•’-*2] 


Along  with  the  definition  of  u^, 


dR 

dt 


-  u 


f  ’ 


(A-24) 


(A-25) 


Equations  (A-9),  (A-ll),  (A-20),  (A-21),  (A-23),  and  (A-25)  are 

solved  to  determine  p,  Ht>  H^,  V,  P^,  ad  Pt  ahen  >  H^. 

The  constants  a  ,  b  ,  d  ,  e  ,  and  e  are  assigned  values  1.3, 
v  v  v  v 

1.0,  0.64,  20.,  and  0.59,  respectively,  based  on  analysis  of  the 
still-air  HTAG  release  experiments  of  Havens  and  Spicer  (Reference 

3). 
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b.  Initial  Gravity  Current  Development 

To  model  the  initial  cloud  shape,  the  tail  and  head  height  are 
considered  constant  with  respect  to  radius.  The  momentum  balance  on 
the  cloud  is  then  given  by 


d_ 

dt 


[  \  *  Pt  ]  -  [  "h  "t  ♦  avbv  ^  ] 


2  dPv 


ir&  d  p  RH.  u_  -  -t— 
v  v'a  n  f  dt 


(A-26) 


where  the  first  term  on  the  right-hand  side  represents  the  static 
pressure  force  on  the  head  and  the  second  term  represents  the  drag 
force  on  the  bottom  surface  of  the  cloud.  The  third  force  is  the 
acceleration  reaction  by  the  ambient  fluid,  represented  by  Equation 
(A-8) . 

The  dimensions  of  the  head  are  again  given  by 

\  =  *  '  avbvHb  (A-27) 

and 


[g&P/p&] 


(A- 28) 


When  the  height  of  the  tail  is  assumed  uniform  with  respect  to 
radius,  it  follows  that 


H 


!  -  ’*l\ (K  *  v 


(A-29) 


where  M  is  the  total  mass  of  the  cloud.  The  momentum  of  the  head  P. 

r 

and  tail  P  are  then 


.  2  ^  <R3  -  ^ 
ph  -  3  - - - uf 


(A-30) 


and 


(A-31) 


A  £  _ 


Pt“  3  * 


PW 


—  uf 


Equations  (A-26)  through  (A-31)  determine  the  momentum  of  the 
blanket  as  a  function  of  time,  and  thus  the  frontal  velocity  u^. 
The  cloud  accelerates  from  rest  because  -  0  initially. 

3.  Material  and  Energy  Balances 

The  balance  on  the  total  mass  of  gas  in  the  source  blanket 
(M  -  wR2H p)  is 


“  -  s  [  *r2h'  ]  -  E<*> +  K  *  K 


(”r2) 


(A-32) 


where  E(t)  is  the  gas  evolution  rate  from  the  primary  (liquid) 

source.  For  spills  over  water,  the  water  entrainment  term  (M  )  is 

w,  s 

included  in  the  source  blanket  description  and  is  calculated  from 
Equation  (A-46) ,  and  the  (humid)  air  entrainment  rate  (Equation 
(A-12))  is 


M  -  2jtRH(£U,)p 
4  I 


•/( 


gApH/(p 


(A-33) 


The  balance  on  the  mass  of  contaminant  in  the  source  blanket 

(M  -  w  *R2Hp)  is 
c  c 


-  3t  [  *Xa'’  ]  - 


(A-34) 


and  the  mass  balance  on  the  air  in  the  source  blanket 


(M  -  w  wR  Hp)  is 

A  B 


d  [  „2„  1  a  ^*max 

-2-  w  jtR  Hp  -  7—---  -  -  w 

dtla  ll+H  w  a 

J  a  c 


(*R2) 


(A-35) 


where  the  ambient  humidity  is  and  the  mass  fraction  of 
contaminant  and  air  are  w  -  M  /M  and  w  -  M  /M,  respectively. 

C  C  B  B 


2 

The  energy  balance  on  Che  source  blanket  (hxR  Hp)  gives 

|  hwR2Hp  1  -  h_E(t)  +  h  M  +  h  H 
dt  L  J  Ev/  aa  ww,s 

.  ^*max  .  _2.  • 

-  h  -  (xR  )  +  Q  (A-36) 

w  s 

c 

where  h  is  the  enthalpy  of  the  emitted  gas,  h  is  the  enthalpy  of 

the  ambient  humid  air,  and  h^  is  the  enthalpy  of  any  water  vapor 

entrained  by  the  blanket  if  over  water.  There  are  three  alternate 

submodels  included  for  the  heat  transfer  (Q  )  from  the  surface  to 

s 

the  cloud. 

The  simplest  method  for  calculating  the  heat  transfer  between 
the  substrate  and  the  gas  cloud  is  to  specify  a  constant  heat 
transfer  coefficient  for  the  heat  transfer  relation 

K  -  %  [  *  (  *2  •  Rp )  ]  -  vT  [  *  (  *2  •  Rp )  ]  <A-37> 

where  Q  is  the  rate  of  heat  transfer  to  the  cloud,  q  is  the  heat 
s  s 

flux ,  and  AT  is  the  temperature  difference.  For  the  calculation  of 
heat  transfer  over  the  source,  the  temperature  difference  is  based 
on  the  average  temperature  of  the  blanket. 

In  the  evaluation  of  the  Burro  and  Coyote  series  of  experiments, 
Koopman  et  al.  (Reference  27)  proposed  the  following  empirical  heat 
transfer  coefficient  relationship  for  heat  transfer  between  a  cold 
LNG  cloud  and  the  ground 

ho  -  Vcp  <A-38) 

where  the  value  of  VH  was  estimated  to  be  0.0125  m/s.  This  constant 
can  be  varied  in  the  model. 

From  the  heat  transfer  coefficient  descriptions  for  heat 
transfer  from  a  flat  plate,  the  following  relationships  can  be 
applied.  For  natural  convection,  the  heat  transfer  coefficient 
is  estimated  using  the  Nusselt  (Nu),  Grashoff  (Gr)  and  Prandtl 
(Pr)  numbers,  McAdams  (Reference  28),  from 


'4’S 


Nu  -  0.14  (Gr  Pr)!/^ 


(A-39) 


2.3 

gP  C  u 


(A-40) 


where  h  is  the  heat  transfer  coefficient  due  to  natural  convection 
n 

and  Pr  is  the  Prandtl  number.  In  order  to  simplify  the 
calculations,  the  parameter  group 


Pr'2  ( V"  ]3  (  h£-t  ) 


(A-41) 


is  estimated  to  be  60  in  mks  units .  The  actual  value  of  the  group 
is  47.25,  58.5,  and  73.4  for  air,  methane,  and  propane,  respec¬ 
tively.  Equation  (A-40)  becomes 

r  -jl/3 

h„-18  [  (fc]24TJ  <A'42) 

where  the  density  p,  molecular  weight  MW,  and  temperature  difference 
AT  are  based  on  the  average  composition  of  the  gas  blanket. 

For  forced  convection,  the  Colburn  analogy  (Reference  29)  is 
applied  to  a  flat  plate  using  the  Stanton  number  for  heat  transfer 
St  and  the  Prandtl  number  as 

H 


StH  Pr 


2/3  ff  J* 

"  2  - 

u 


(A-43) 


hf  -  (upCp)  Pr 


-2/3  U* 


(A-44) 


-■  -  A'V-V  ■ V'.V- 


where  is  the  heat  transfer  coefficient  due  to  forced  convection. 
If  the  velocity  is  evaluated  at  the  top  of  the  gas  blanket  and  Pr  is 
estimated  to  be  0.741, 


u*  zn 

1.22  —  ^ 

"o  H 


(*-45) 


The  overall  heat  transfer  coefficient  is  then  the  maximum  of  the 

forced  and  natural  coeff icients ,  i.e.  hrt  -  max(h_,h  ).  The  heat 

0  f  n 

flux  and  transfer  rate  are  then  estimated  by  Equation  (A-37). 

If  the  gas  blanket  is  formed  over  water,  water  will  be 
transferred  from  the  surface  to  the  cloud  by  a  partial  pressure 
driving  force  associated  with  the  temperature  difference  between  the 
surface  and  the  gas  blanket.  The  rate  of  mass  transfer  of  water  is 


( C  -  p».c )  [  *  (  *2  •  Rp  ]  ] 


(A-46) 


where  Fq  is  the  overall  mass  transfer  coefficient.  The  driving 
force  is  the  difference  of  the  vapor  pressure  of  water  at  the 
surface  temperature  p  and  the  vapor  pressure  of  water  at  the 


cloud  temperature  p 


The  natural  convection  coefficient  is  based 


on  the  heat  transfer  coefficient  and  the  analogy  between  the 
Sherwood  number  (Sh)  and  the  Nusselt  number  (Nu)  suggested  by  Bird 
et  al.  (Reference  30) 


Sh  -  0.14  (Gr  Sc) 


FnL  f  WJ  1 

*  *  I  P  J 


(A-47) 


If  the  Schmidt  number  is  taken  as  0.6,  and 


M-l 

(t  mwJ 


is  estimated  to 


be  2.2  x  10  in  mks  units, 


F  -  9.9  x  10 
n 


(M: 


(A-48) 


V  -•  1 


For  forced  convection,  Treybal  (Reference  29)  suggests  that  the 
Stanton  number  for  mass  transfer  StM  and  the  Stanton  number  for  heat 

n 

transfer  Stu  are  related  by 
n 

StM  *  StH  (  H  ]V3-  115  StH 
Or, 


20.7  h. 


MW  C 
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(A-50) 


The  overall  mass  transfer  coefficient  Fq  is  calculated  as  the  larger 

of  the  natural  and  forced  convection  coefficients. 

For  the  case  when  the  primary  (liquid)  source  emission  rate  E(t) 

2 

is  larger  than  the  atmospheric  takeup  rate  Equations  (A- 

32),  (A- 34) ,  (A-35),  and  (A-36)  are  integrated  for  the  mass, 

concentration,  and  enthalpy  of  the  gas  blanket  along  with  an 

appropriate  equation  of  state  (i.e.r relationship  between  enthalpy 

and  temperature  and  between  temperature  and  density). 

For  the  case  when  the  emission  rate  is  not  sufficient  to  form  a 

gas  blanket,  the  flux  of  contaminant  is  not  determined  by  the 

maximum  atmospheric  takeup  rate.  Consider  the  boundary  layer  formed 

by  the  emission  of  gas  into  the  atmosphere  above  the  primary  source. 

If  the  source  is  modeled  to  have  a  uniform  width  2b  and  entrain  no 

air  along  the  sides  of  the  layer,  the  balance  on  the  total  material 

(p  u  H  )  in  a  differential  slice  of  the  layer  is 
L  L  L 

-t-  P-U.H,  •»  p  w  + 
dx  [  L  L  L  J  *a  e 


(A-51) 


where  w  is  the  vertical  rate  of  air  entrainment  into  the  layer 
e 

given  by  Equation  (A-83),  p^  is  the  average  density  of  the  slice, 

and  Q*/»c  is  the  total  flux  of  gas  from  the  primary  (liquid)  source. 

The  balance  on  the  mass  flow  rate  of  contaminant  (w  p  u  H  )  at  any 

C  L  L  L 

(X  -  X  )  is 
up' 


c  u  H.  - 
C  ,  L  Lt  L 


(A-52) 
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With  an  equation  of  state  to  relate  cc  ^  and  p^,  Equation  (A-51)  is 

integrated  from  the  upwind  edge  of  the  source  (x  -  x  )  t0  t*ie 

downwind  edge  (x  -  L  +  x  ) . 

up 

In  order  to  generate  the  initial  conditions  for  the  downwind 
dispersion  calculations,  the  maximum  concentration  cc  and  the 
vertical  dispersion  parameter  S  are  needed.  Since  Equations  (A-51) 
and  (A-52)  are  written  for  a  vertically  averaged  layer,  consider  the 
vertical  average  of  the  power  law  distribution.  The  height  of  the 
layer  H.  is  the  height  to  some  concentration  level,  say  10  percent  of  the 
maximum.  Although  strictly  a  function  of  a,  this  value  is  modeled 
by 


HL  “  5l  heff 


(A-53) 


where  HgFF  is  the  effective  height  defined  by  Equation  (A-79)  and  6^ 
is  2.15.  The  vertically  averaged  concentration  c  _  can  be  defined 

CtL 

by 


:.L  "L  -  £** 


(A- 54) 


And  similarly,  the  effective  transport  velocity  u^  is  defined  by 

f cu  dz  (A- 55) 


'c.lVl  J  q  x 


With  Equation  (A-53)  and  defining  relations  for  H£FF  and  uEFp 
(Equations  (A-79)  and  (A-93),  respectively),  it  follows  that 

(A-56) 


and 
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(A-57) 


where  w'  is  given  by  Equation  (A-83). 
€ 


(A-58) 


4.  Maximum  Atmospheric  Takeup  Rate 


The  maximum  atmospheric  takeup  rate  will  be  the  largest  takeup 
rate  which  satisfies  Equations  (A-51)  and  (A-52).  As  well,  the 
maximum  concentration  of  contaminant  In  the  power  law  profile  at  the 
downwind  edge  of  the  source  will  be  the  source  contaminant 
concentration  (c  )  .  If  Equations  (A-51)  and  (A-52)  are  combined 
along  with  the  assumption  of  adiabatic  mixing  of  Ideal  gases  with 
the  same 

P  -  Pa 

constant  molal  heat  capacity  (l.e.  -  -  7  -  constant),  the 

c 

maximum  takeup  flux  Is  modeled  by 


(A-59) 


(A-60) 


An  upper  bound  of  the  atmospheric  takeup  flux  can  be 
characterized  by  the  condition  where  the  source  begins  to  spread  as 
a  gravity  Intrusion  against  the  approach  flow.  In  water  flume 
experiments,  Brltter  (Reference  22)  measured  the  upstream  and 
lateral  extent  of  a  steady  state  plume  from  a  circular  source  as  a 
function  of  Ri#.  A  significant  upstream  spread  was  obtained  for  Rl^ 
>  32,  and  lateral  spreading  at  the  center  of  the  source  was 
Insignificant  for  Ri^  <  8.  The  presence  of  any  significant  lateral 
spreading  represents  a  lower  bound  on  the  conditions  of  the  maximum 
takeup  flux. 

The  integral  of  Equation  (A-60)  is  calculated  using  a  local 
Richardson  number  of 


Ri*(x)  - 


f(x  -  x  ) 
*  up' 
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where 
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(A-62) 

and  is  3.1  (corresponding  to  Ri^  -  20(8  <  Ri^  <  32)).  Using  this 
Ri^(x)  and  the  first  two  terms  of  ^(Ri^),  Equation  (A-60)  is 
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dx 


1.04 

-  ..  .  .1.04  l+o 

0.88  +  0.099  f  x 


In  order  to  simplify  the  numerical  problem,  the  integral  is 
approximated  as 


1 
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.1.04 


In 


0.88  +  0.099  f 


1.04  . l l+o  J 


0.88 


<f>  0.099  L  f 

which  then  specifies  the  maximum  atmospheric  takeup  flux. 

5.  Transient  HTAG  Release  Simulation 


(A-63) 


If  a  steady  state  spill  is  being  simulated,  the  transient 
source  calculation  is  carried  out  until  the  source  characteristics 
are  no  longer  varying  significantly  with  time.  The  maximum 
centerline  concentration  cc>  the  horizontal  and  vertical  dispersion 
parameters  and  S^,  the  half-width^ b,  and  if  necessary,  the 
enthalpy,  h, are  used  as  initial  conditions  for  the  downwind 
calculation  specified  in  a  transient  spill. 

If  a  transient  spill  is  being  simulated,  the  spill  is  modeled 
as  a  series  of  pseudo- steady-state  releases.  Consider  a  series  of 
observers  traveling  with  the  wind  over  the  transient  gas  source 
described  above;  each  observer  originates  from  the  point  which 
corresponds  with  the  maximum  upwind  extent  of  the  gas  blanket 
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(x  -  -R  ) .  The  desired  observer  velocity  is  the  average  transport 
max 

velocity  of  the  gas  u^p  from  Equation  (A- 93);  however,  the  value  of 

u  _  will  differ  from  observer  to  observer  with  the  consequence  that 
£FF 

some  observers  may  be  overtaken  by  others.  For  a  neutrally  buoyant 

cloud,  u___  becomes  a  function  of  downwind  distance  alone  which 
EFf 

circumvents  this  problem.  With  this  functionality,  Colenbrander 
(Reference  17)  models  the  observer  velocity  as 


ut(x) 


1  +  a  I 


,a/(l+a) 


x  +  R 


R  +  R 
2  m  max 


(A-64) 


where  is  the  value  of  S£  when  the  averaged  source  rate 

(*R2Q#)  ?s  a  maximum  and  the  subscript  i  denotes  observer  i.  Noting 
that  u^(x)  -  dXj/dt,  observer  position  and  velocity  as  functions  of 
time  are  determined. 

A  pseudo- steady estate  approximation  of  the  transient  source  is 
obtained  as  each  observer  passes  over  the  source.  If  t  and  t. 


denote  the  times  when  Observer  i  encounters  the  upwind  and  downwind 
edges  of  the  source  respectively,  then  the  source  fetch  seen  by 
Observer  i  is: 


L.  -  x  -  x. 

1  Upi  dni 


(A-65) 


The  width  of  the  source  2B^(t)  is  defined  by 


B'2(t)  -  R2(t)  -  xj(t) 


(A-66) 


Then  the  gas  source  area  seen  by  Observer  i  is 


r  dnt 

V>i-2  Jt  W 


(A-67) 


where  2b^  is  the  average  width. 


The  takeup  rate  of  contaminant  (Q*Lb)^  is  calculated  as 
tdn 

2<v*>i  - 2  {  1  viv* 

tui>i 

The  total  mass  flux  rate  from  the  source  is 


o  ftdni 

V 

2<'LULHLb>i  -  2  J, 

p  w'  + 
a  e 

w 

i  °  j 

4 

(A-68) 


(A-69) 


With  these  equations,  the  average  composition  of  the  layer  can 
be  determined  at  each  x 
the  layer  given  by 


x  over  the  source.  With  the  enthalpy  of 


(hT 


Q*’ 

1  h 

V 

J  t 

iin 

l  c 

BiUidt 


(A-70) 


due  to  the  choice  of  the  reference  temperature  as  the  ambient 

temperature  along  with  a  suitable  equation  of  state  relating 

enthalpy,  temperature,  and  density,  the  source  can  be  averaged  for 

each  observer.  After  the  average  composition  of  the  layer  is 

determined  at  the  downwind  edge,  an  adiabatic  mixing  calculation  is 

performed  between  this  gas  and  the  ambient  air  when  applicable. 

This  calculation  represents  the  function  between  density  and 

concentration  for  the  remainder  of  the  calculation  if  the 

calculation  is  adiabatic;  it  represents  the  adiabatic  mixing 

condition  if  heat  transfer  is  included  in  the  downwind  calculation. 

For  each  of  several  observers  released  successively  from 

x  =*  -R  ,  the  observed  dimensions  L  and  b,  the  downwind  edge  of  the 
max 

source  x^,  the  average  vertical  dispersion  coefficient  S^,  the 

average  takeup  flux  Q. ,  the  centerline  concentration  c  ,  and  if 

w  c 

applicable,  the  average  enthalpy  h  can  be  determined  for  each 

la 

observer.  With  these  input  values,  a  steady  state  calculation  is 
made  for  each  observer.  The  distribution  parameters  for  any 
specified  time  ts  are  determined  by  locating  the  position  of  the 
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series  of  observers  at  time  t  ,  i.e.  x.(t  ).  The  corresponding 

S  1  s 

concentration  distribution  is  then  computed  from  the  assumed 
profiles . 


B.  STEADY-STATE  DOWNWIND  DISPERSION 

The  model  treats  dispersion  of  gas  entrained  into  the  wind 

field  from  an  idealized  rectangularly  shaped  source  of  width  2b  and 

length  L.  The  circular  source  cloud  is  represented  as  an  equivalent 
2  2 

area  square  (L  -  itK  -  2bL) .  Similarity  forms  for  the 
concentration  profiles  are  assumed  which  represent  the  plume  as 
being  composed  of  a  horizontally  homogeneous  section  in  which  only 
vertical  dispersion  is  considered,  with  Gaussian  concentration 
profile  edges  as  follows: 

c(x,y,z)  -  c  (x)  exp 
c 


-  cc<x)  exp 


|y|  -  Mx) 

2  f  2  1 

1+a 

S  (x) 

y  J 

l  z  J 

for  | y |  >  b 
for  | y |  <  b 


S2(x) 


1+Ot 


A  power  law  wind  velocity  profile  is  assumed 


(A-71) 


(A-72) 


where  the  value  of  a  is  determined  by  a  weighted  least  squares  fit 
of  the  logarithmic  profile 


(A-73) 


Functional  forms  for  #  and  typical  values  of  a  are  given  In  Table 
A-l  for  different  Pasquill  stability  categories.  With  these 
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profiles,  the  parameters  of  Equation  (A- 71)  are  constrained  by 
ordinary  differential  equations. 

1.  Vertical  Dispersion 


The  vertical  dispersion  parameter  is  determined  by  requiring 


that  it  satisfy  the  diffusion  equation 


3c  d  v  dc 
u  —  *  —  K  — 
x  3x  3z  z  3z 


with  the  vertical  turbulent  diffusivity  given  by 


(A-74) 


Kz  *(RiJ 


(A-75) 


The  function  <£(Ri#)  is  a  curve  fit  of  laboratory  scale  data  for 
vertical  mixing  in  stably-density-stratified  fluid  flows  reported  by 
McQuaid  (Reference  21),  Kantha  et  al.  (Reference  29),  and  Lofquist 
(Reference  20)  for  Ri^  >  0.  For  Ri^  <  0,  the  function  <£(Ri^,)  is 
taken  from  Colenbrander  and  Puttock  (Reference  18)  and  has  been 
modified  so  the  passive  limit  of  the  two  functions  agree  as  follows: 


4(Ri*)  -  0.88  +  0.099  Ri*'04  +  1.4  x  10'25  Ri*’7 


-  0. 88/(1  +  0.65  { Ri . |  ) 


Ri.  >  0 
* 


Ri.  <  0 
* 


(A-76) 


The  friction  velocity  is  calculated  using  Equation  (A-73)  from 
a  known  velocity  u^  at  a  specific  height  Zq.  Combining  the  assumed 
similarity  forms  for  concentration  and  velocity,  Equations  (A-71), 
(A-72) ,  (A-74),  and  (A-75)  give 


ku*(l  +  o) 

*<RiJ 


where  the  Richardson  number  Ri*  is  computed  as 


(A-77) 


Ri*  -  8 


p  '  Pa  1  HEFF 


(A-78) 


and  the  effective  cloud  depth  is  defined  as 


EFF  c 


i-  fed*  .  rf  1  1  Ji- 
ccJ0  U  +  «  J  1  +  ° 


(A- 79) 


Equation  (A- 77)  can  be  viewed  as  a  volumetric  balance  on  a 
differential  slice  of  material  downwind  of  the  source.  For  a  mass 


balance  over  the  same  slice. 


[  VA  ]  -  V-e 


(A-80) 


which  is  the  same  result  as  Equation  (A-51)  without  the  source  term. 
With  Equations  (A-57)  and  (A-58),  this  is 


dx  (  PLUEFFHEFF  )  paWe 


(A-81) 


Using  the  assumption  of  adiabatic  mixing  of  ideal  gases  with 


the  same  constant  molal  heat  capacity  (i.e. 


P  '  P « 


constant)  along 


with  the  contaminant  material  balance,  the  mass  balance  becomes 


dx  (  UEFFHEFF  )  ”  W€ 


(A-82) 


which  leads  to 


v  hn  (1  +  q) 


e  6, 


*(Ri*) 


(A-83) 


Equations  <A-81)  and  (A-83)  are  combined  to  give 


(  PLUEFFHEFF  ) 


^ku^l  +  a) 
*(Ri*) 


(A-84) 


Furthermore,  Equation  (A-84)  is  assumed  to  apply  when  (p  -  p  )/c 

&  c 

is  not  constant. 


When  heat  transfer  from  the  surface  is  present,  vertical  mixing 
will  be  enhanced  by  the  convective  turbulence  due  to  heat  transfer. 
Zeman  and  Tennekes  (Reference  33)  model  the  resulting  vertical 
turbulence  velocity  as 


w  .  1  * 

—  **  1  +  7  — 

u  4  u . 

*  * 


where  w^  is  the  convective  scale  velocity  described  as 


(A-85) 


v  T  -  T  _ 
gh  l  s  c ,  L 

-  T  , 


(A-86) 


If  u  is  evaluated  at  H  , 

Lb  r 


2-  -  r  1 4  i 

u.  L  * 


(A-87) 


where 


RiT  “  8 


Ts  "  Tc,L  1  HEFF 


v*0  heff 


(A-88) 


and  T  is  the  temperature  obtained  from  the  energy  balance  of 
C ,  L* 

Equations  (A-102)  and  (A-103).  Equation  (A-84)  is  modified  to 
account  for  this  enhanced  mixing  by 


p  kw(l  +  a) 


d  f  _  1  ^aK  : 

dx  [  PL  EFF  EFF  J  *(ri») 


(A-89) 


where  Ri'  -  Ri .  — 

*  *  w 


Although  derived  for  two-dimensional  dispersion,  this  is 
extended  for  application  to  a  HTAG  plume  which  spreads  laterally  as 
a  density  intrusion: 


with  the  horizontal  turbulent  diffusivity  given  by 


Ky  “  K0UxBEFF 

When  b  -  0,  S  -  Jl  a  ,  where  a  is  the  similarity  parameter 

y  y  y  a 

correlated  by  Fasquill  (Reference  31)  in  the  form  a -  6x  . 
Furthermore,  Equations  (A-95)  and  (A-96)  require  that 

d  a  y, 

„  __y  .  v 

y  dx  0  EFF 


(A-96) 


(A-97) 


where  7-  -  2  -  1/0  and  Krt  - 

L  U  X 

rl/^ 


Then 


S  ^-^B2 
y  dx  ir  EFF 


beff 


(A- 98) 


where  Equation  (A-98)  is  also  assumed  applicable  for  determining 
when  b  is  not  zero. 

At  the  downwind  distance  xfc  where  b  «  0,  the  crosswind 
concentration  profile  is  assumed  Gaussian  with  given  by 


S  -  Jl  6(x  +  x  ) 
y  P 


P 


(A-99) 


vhere  x  is  a  virtual  source  distance  determined  as 
v 


S  (x  )  -  Jl  6(x  +  x  )/9 

y  t  v  t  p 


(A-100) 


The  gravity  spreading  calculation  is  terminated  for  x  >  x  . 

For  a  steady  plume,  the  centerline  concentration  cc  is 
determined  from  the  material  balance 


U-Z. 

r  s 

1  cu  dydz  -  2c 

•»o  J-«  x  c 

0  0 

z 
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(A- 101) 


where  E  is  the  plume  source  strength. 


3 .  Energy  Balance 


For  some  simulations  of  cryogenic  gas  releases,  heat  transfer 
to  the  plume  In  the  downwind  dispersion  calculation  may  be 
Important,  particularly  In  low  wind  conditions.  The  source 
calculation  determines  a  gas/air  mixture  initial  condition  for  the 
downwind  dispersion  problem.  Air  entrained  into  the  plume  is 
assumed  to  mix  adiabatically .  Heat  transfer  to  the  plume  downwind 
of  the  source  adds  additional  heat.  This  added  heat  per  unit  mass 
is  determined  by  an  energy  balance  on  a  uniform  cross  section  as 

dx  [  DhPLUEFFHEFF  ]  "  qs^SL  (A-102) 

where  q  is  determined  by  Equation  (A-37)  along  with  the  desired 
s 

method  of  calculating  h^.  Equation  (A-102)  is  applied  when  b  «■  0 
and  is  extended  to 


dx  [  DhPLUEFFHEFFBEFF  ]  “  qs  ®EFF  '  6L  (A- 103) 

when  b  >  0.  Since  the  average  density  of  the  layer  p L  cannot  be 
determined  until  the  temperature  (i.e.  D^)  is  known,  a  trial  and 
error  procedure  is  required. 

Equations  (A-77),  (A-78),  (A-79),  (A-87)-(A-91) ,  (A-94) ,  and 

(A-98)-(A~103)  are  combined  with  an  equation  of  state  relating  cloud 

density  to  gas  concentration  and  temperature  and  are  solved 

simultaneously  to  predict  S  ,  S  ,  c  ,  and  b  as  functions  of  downwind 

z  y  c 

distance  beginning  at  the  downwind  edge  of  the  gas  source. 


C.  CORRECTION  FOR  ALONG-WIND  DISPERSION 


Following  Colenbrander  (Reference  17),  an  adjustment  to  the 

values  of  cc  is  applied  to  account  for  dispersion  parallel  to  the 

wind  direction.  The  calculated  centerline  concentration  c  (x)  is 

c 

considered  to  have  resulted  from  the  release  of  successive  planar 
puffs  of  gas  (cc(x)Ax)  without  any  dispersion  in  the  x-direction. 
If  it  is  assumed  that  each  puff  diffuses  in  the  x-direction  as  the 


puff  moves  downwind  independently  of  any  other  puff  and  that  the 
dispersion  is  one-dimensional  and  Gaussian,  the  x-direction 
concentration  dependence  is  given  by 


■s 
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CM 


c :  )  - 


c  <x  )Ax 
C  Pj  i 


(A-104) 


where  x  denotes  the  position  of  the  puff  center  due  to  observer  i. 


After  Beals  (Reference  34),  the  x-direction  dispersion 


coefficient  is  assumed  to  be  a  function  of  distance  from  the 


downwind  edge  of  the  gas  source  (X  -  x  -  x^)  and  atmospheric 
stability  given  by 


ax(X)  -  0.02  X* 


unstable,  x  >  130  m 


0.0 4  X 


neutral,  x  >  100  m 


(A-105) 


-  0.17  Xv 


stable,  x  >  50  m 


where  (X  •  x  -  x^)  and  are  in  meters.  The  concentration  at  x  is 


then  determined  by  superposition,  i.e.,  the  contribution  to  c  at  a 


given  x  from  neighboring  puffs  is  added  to  given  an  x-direction 


corrected  value  of  c^.  For  N  observers. 


Cc<x)  “ 


M*  ) 


6XP  -  2 


(A- 106) 


and  for  large  N, 


i  r  cc(° 
'C°°  "  JR  Jo  V*  -  V 


1  x  -  i 


2  *x<€  -  V 


(A-107) 


The  corrected  centerline  concentration  c£  is  used  in  the 
assumed  profiles  in  place  of  cc>  along  with  the  distribution 
parameters  Sy,  Sz,  and  b. 

D.  DEGADIS  MODEL  INPUTS  AND  OUTPUTS 

As  implemented  under  VAX/VMS?  DEGADIS  requires  three  areas 
of  input  information: 

VAX/VMS  command  procedure  for  execution 
simulation  definition 
numerical  parameters 

The  VAX/VMS  command  procedure  used  to  execute  DEGADIS  is 
generated  in  DEGADISIN  by  default.  As  well,  DEGADISIN  is  the 
interactive  input  module  which  generates  the  simulation 
definition  from  user  responses.  An  example  input  session  is 
included  in  Section  3.  The  numerical  parameters  (convergence 
criteria,  initial  increments,  etc.)  are  supplied  to  DEGADIS 
through  a  series  of  input  files.  Although  these  numerical 
parameters  are  easily  changed,  the  user  should  need  to  change 
these  only  rarely  with  the  exception  of  the  time  sort 
parameters.  Additional  information  can  be  found  in  Reference  3. 

1.  VAX/VMS  Command  Procedure 

The  VAX/VMS  command  procedure  generated  by  DEGADISIN 
controls  the  execution  of  images  for  the  simulation.  Image 
execution  follows  one  of  two  paths,  either  for  a  transient 
release  or  for  a  steady  state  release.  DEGADISIN  will 
automatically  generate  the  appropriate  command  procedure. 
DEGADISIN  requires  a  simulation  name  be  specified.  The 
simulation  name  must  be  a  valid  VAX/VMS  file  name  without  a  file 
extension  and  is  designated  herein  as  RUNNAME.  DEGADIS  will  use 

*VAX  and  VMS  are  registered  trademarks  of  Digital  Equipment 
Corporation. 


this  file  name  with  standard  extensions  for  input,  interprocess 
communication,  and  output. 

2.  Simulation  Definition 

DEGADISIN  is  an  interactive  method  of  simulation  definition 
where  the  user  specifies  information  about  the  ambient  wind 
field,  the  properties  of  the  released  gas,  and  some  details  of 
the  release.  A  summary  of  required  input  information  is 
included  in  Figure  A- 5. 

The  ambient  wind  field  is  characterized  by  a  known  velocity 
Uq  at  a  given  height  zq,  a  surface  roughness  zR,  and  the 
Pasquill  stability  class.  The  Pasquill  stability  class  is  used 
to  estimate  values  of  the  lateral  similarity  parameter 
coefficients  6  and  p  (Reference  31),  values  of  the  along-wind 
similarity  coefficients  (Reference  34) ,  and  the  Monin-Obukhov 
length  A  used  by  Businger  et  al.  (Reference  32)  in  their 
logarithmic  velocity  profile  function.  The  Monin-Obukhov  length 
is  then  used  to  calculate  the  friction  velocity  u*.  Once  these 
parameters  have  been  estimated  using  the  Pasquill  stability 
class,  the  user  has  the  option  of  interactively  changing  any  of 
these  to  better  describe  the  simulation.  In  addition  to  these 
specifications,  the  ambient  temperature,  pressure,  and  humidity 
must  be  specified. 

The  properties  of  air  and  the  released  gas  are  used  to 
evaluate  the  mixture  density  as  a  function  of  temperature  and 
composition.  The  desired  released  gas  properties  include  the 
molecular  weight  MWC,  the  storage  temperature  (normal  boiling 
point  for  cryogenic  gases)  Tg,  the  vapor  phase  density  at  the 
storage  temperature  and  ambient  pressure  pq,  and  two  constants 
q^  and  p^  which  describe  the  heat  capacity  according  to  the 
equation 
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Figure  A-5.  Summary  of  simulation  definition  input 
information  for  DEGADIS. 
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Figure  A-5.  (continued) 
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Variable  Symbol 

CAS_crr 
GAS_UFL 
CAS_m 
GAS  1ST 


n/a  tan  for  coot— In— t  b— t 

capacity 

aol*  Upper  (In— ablllt?  limit 

fraction 

molt  lowr  flaooa billty  limit 

fraction 

n  Height  for  coo toura 


Variable 

ar 

snu.D 

HHtf.I) 

m(3.i) 


ISOTLdO 

lootbornml  Ssleae 
or  btiml 
Adiabatic  Mixing 
Calculation 


Number  of  denalty  points 
Contaminant  mole  fraction 

Contaminant  concentration 
Mixture  desalt? 


cncM-.rALSs. 


Stood?  State  Sol— 


k* )• 

Stood?  etate  reloaaa  —to 

„  ■ 

Stead?  etate  p rimer?  source 

P 

radlua 

Figure  A-5.  (continued) 
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C  (T)  -  (MW,)"1  3.33  x  104  +  qx 


P1  P1 
T  L  -  T 

1  XQ 
T  -  T_ 


(A-108) 


where  C  (T)  is  the  mean  heat  capacity  (J/kg  K)  at  temperature 
pc 

T.  Note  that  a  constant  heat  capacity  with  respect  to 
temperature  can  be  obtained  by  setting  —  1.0  and  choosing 
the  appropriate  value  for  q^.  Representative  gas  properties 
for  liquefied  natural  gas  (LUG)  as  methane,  liquefied 
petroleum  gas  (LPG)  as  propane,  pure  unreacted  NO2 ,  and  pure 
unreacted  N2O4  are  included  in  DEGADISIN.  Also  included  are 
the  lower  and  upper  flammability  limits  (LFL  and  UFL, 
respectively)  for  LNG  and  LPG. 

The  user  may  also  choose  to  calculate  the  mixture 
density  as  a  function  of  composition  using  some  other 
method.  This  mixture  density  is  entered  in  the  program  as 
if  the  release  were  isothermal;  for  each  composition,  the 
program  requests  the  contaminant  mole  fraction,  the 
contaminant  concentration,  and  the  mixture  density.  For 
ease  of  input,  these  values  may  be  entered  from  a  file  made 
available  to  DEGADISIN. 

In  specifying  the  details  of  the  release,  the  user  must 
choose  to  simulate  the  release  as  transient  or  steady  state. 
For  both  release  types,  the  area  source  is  assumed  circular. 
The  source  radius  and  emission  rate  must  be  specified  for  a 
steady  state  release  only  once,  while  these  must  be 
specified  as  a  function  of  time  for  transient  releases 
(either  interactively  or  by  file).  For  transient  releases, 
the  user  must  specify  the  initial  amount  of  gas  present  over 
the  source  (in  order  to  simulate,  for  example,  instantaneous 
releases  such  as  the  Thorney  Island  Trials). 

Figure  A-5  summarizes  the  simulation  information 
gathered  by  DEGADISIN  contained  in  the  RUNNAME  file  with 
extension  INP.  The  structure  of  RUNNAME. INP  is  illustrated 
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in  Figure  A- 6.  At  this  point,  RUNNAME . INP  may  be  edited  to 
correct  any  misinformation  entered  during  the  input  session. 
Note  that  care  must  be  exercised  when  editing  RUNNAME. INP 
due  to  the  fact  that  information  contained  in  the  file  can 
be  different  depending  on  the  answered  questions  (e.g. 
steady  state  versus  transient  simulation). 


EAGLE  6  test  situation  as  an  ISOTHERMAL  release  of  PURE  M02  WITHOUT 
CHEMICAL  REACTION  usinS  the  HIGH  SOURCE  RATE.  The  source  Sas  is  assuaed  to 
have  the  saae  temperature  as  the  LIQUID  POOL  TEMPERATURE  (15C). 

5.580000  12.00000  1 . 0C00000E-06 

4 

0.1300000  0.9000000  O.OOOOOOOE+OO 

4 .000C000E-02  1.140000  100.0000 

295.7500  0.9090000  6.7578776E-03 

1  295.7500 

0  0.0000000E+00 
0  O.OOOOOOOE+OO 
«J2 

46.00000  295.7500  1.769000 

3845,000  1.000000 

1 .00000 OOE-03  5.0000000E-04  0.5000000 

o 

O.0000D00E+00  0.0000000E+00  1.081600  O.OOOOOOOE+OO  295.7500 

1.000000  1.769000  1,769000  O.OOOOOOOE+OO  295.7500 

3  .OOOOO'OL-O-t 
O.OOOOOOOE+OO 
4 

O.OOOCOOOE+OO  1.720000  10.00000 

6023.000  1.720000  10.00000 

5024.000  O.OOOOOOOE+OO  O.OOOOOOOE+OO 
6025.000  O.OOOOOOOE+OO  0.0000000E+00 
F  +  r  F  T  F 

9-fW-l986  14111:56.90 
1.720000  17.  ■’2454  17.72454 


Figure  A-6  .  EAGLE 6 . INP  LISTING . 


3.  Example  Input  Sessions 


The  input  procedures  for  simulation  of  the  transient 
release  and  the  steady  state  release  are  very  similar. 

Therefore,  only  the  specification  of  a  steady- state  release 
(Eagle  3)  has  been  included.  In  the  point-by-point 
discussion  of  the  input  procedure,  note  the  following: 

(*)  The  file  name  specification  RUNNAME  must  satisfy  system 
restrictions . 

(*)  A  line  terminator  (normally  a  carriage  return)  must  end 
every  line  entered  by  the  user. 

(*)  When  DEGADISIN  requests  the  user  to  choose  an  option, 
all  acceptable  responses  are  a  single  character 
(capital  or  lower  case) .  The  default  responses  are 
denoted  by  a  capital  letter  inside  angle  brackets  (e.g. 
<N>) .  When  applicable,  a  menu  of  acceptable  responses 
is  included  inside  parentheses. 

(*)  For  numerical  responses,  a  comma,  space,  tab,  or  line 
terminator  (carriage  return)  may  separate  the  numbers. 

(*)  When  a  file  is  used  as  input  (i.e.,  for  the  density  or 
transient  source  input) ,  DEGADISIN  reads  the  same 
information  from  the  file  which  would  be  entered  at  the 
terminal  in  the  same  order  and  in  the  same  format. 


'TT'JT 


1TI  ’ 


i  RUN  DEGADISIN 


DEnse  GAs  DISpersion  Model  input  aodule. 

Enter  the  siaulation  naae  I  [DIR3RUNNAME  EAGLE6 
INPUT  MODULE  -  DEGADIS  MODEL 


mmmtmtmmmtmmmt 

Enter  Title  Block  —  up  tc  4  lines  of  80  characters 
To  stopi  type  *//' 

EAGLE  6  test  siaulation  as  an  ISOTHERMAL  release  of  PURE  N02  WITHOUT 
CHEMICAL  REACTION  usinS  the  HIGH  SOURCE  RATE.  The  source  das  is  assured  to 
have  the  sane  teaperature  as  the  LIQUID  POOL  TEMPERATURE  (150. 


ENTER  WIND  PARAMETERS  --  UO  (a/s!>  ZO  (»)>  and  ZR(») 
UO  --  Wind  velocity  at  reference  heisht  ZO 
ZR  —  Surface  Roudhness 
5.53fl2.00il.E-6 

Enter  the  Pasouill  stability  class.*  (ArF tCtDtEiF)  <D> 


The  values  for  the  ataospheric  Parareters  are  set  as  follows! 

delta:  o,i3oo 

beta;  0.9000 


Monin-OboRhov  lensth!  infinite 

Siflaa  X  Coefficient:  0,0400 

SiSaa  X  Power!  1.1400 

Sidae  X  Miniaua  Distance!  100,0000  a 

Do  you  wish  to  chanSe  any  of  these? 

(No>Delt.'8et3iLendthiCoefficient»Power»Mini8ua)  <N> 

Enter  the  anbie.it  te»?erature(C)  and  pressure(ata)!  22, 6i  .909 

The  aobient  huaidits  can  be  entered  as  Relative  or  Absolute. 
Enter  cither  R  or  A  <R  or  9>l 
Enter  the  relative  huaidits  <Z)!  35. 

Aabient  Air  density  is  1.0816  KS/a«3 


1 


I 

I 

NOTES  ON  STEADY- STATE  SIMULATION  OF  EAGLE6  ! 


Begin  the  input  procedure  by  execution  of  DEGADISIN. 


The  file  name  specification  must  follow  system 
restrictions.  The  DEGADIS  model  uses  this  file  name 
along  with  various  file  extensions  for  input  and 
output. 

The  Title  Block  is  used  to  carry  any  desired  comments 
such  as  information  on  the  specification  of  certain 
parameters . 


The  wind  field  parameters  include  the  wind  velocity 
(m/s)  at  a  specified  height  (m)  and  the  surface 
roughness  (m) . 


© 

® 


The  Pasquill  stability  class  is  used  to  generate 
estimates  of  other  atmospheric  parameters  which  follow. 

The  current  settings  of  pertinent  atmospheric 
parameters  are  displayed  in  this  list.  If  any  of  these 
are  to  be  changed,  the  first  letter  of  the  parameter  to 
be  changed  is  entered.  Note  that  the  default- - 
indicated  by  <N>--is  No  for  no  changes. 

The  value  of  the  Monin-Obukhov  length  is  indicated  to 
be  infinite  (which  is  the  default  value  for  a  Pasquill 
D  stability).  Any  value  of  the  Monin-Obukhov  length 
can  be  specified,  but  an  infinite  value  is  a  special 
case  which  can  be  obtained  for  any  Pasquill  stability 
class  by  setting  the  Monin-Obukhov  length  to  0. 

The  ambient  temperature  and  pressure  are  entered. 


DEGADISIN  calculates  the  ambient  air  density  for  the 
given  input  parameters. 
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10)  Is  this  an  Isothermal  spiII?  <y  or  M>  Y 


Enter  the  code  name  of  the  diffusing  species!  N02 


The  characteristics  for  the  das  are  set  as  follows! 


Molecular  weight!  46*00 

Storage  temperature  CM!  295.75 

Densits  at  storage  temperature*  PAMB  [kg/m*l33!  1.7231 

Mean  Heat  capacity  constant  3845.0 

Mean  Heat  capacity  power  1.0000 

Upper  Flammability  Limit  [mole  frac3  1.00000E-03 

Lower  Flammability  Limit  [mole  frac3  5.00000E-04 

Height  of  Flammability  Limit  [ml  0.50000 

Dr  you  wish  to  change  any  of  these?  (No*MoleiTe*p*Den*Heat*Power*UpperiLower*Z)  <N>  S 
Enter  the  desired  Density  at  Storage  Temperature  and  ambient  pressure!  1.769 


The  characteristics  for  the  gas  are  set  as  follows! 


Molecular  weight! 

Storage  temperature  [K3! 

Density  at  storage  temperature*  PAHB  Ckg/m*t33! 

Mean  Heat  capacity  constant 

Mean  Heat  capacity  power 

Upper  Flammability  Limit  [mole  frac3 

Lower  Flammability  Limit  [mole  frac! 

Height  of  Flammability  Limit  [m3 


46.00 

295.75 

1.7690 

3845.0 

1.0000 

1.00000E-03 

5.00000E-04 

0.50000 


Do  you  wish  to  change  any  of  these?  (No*HoleiTemp*Den*Heat*Power*UpperrLower*Z)  <H> 


The  density  is  determined  as  a  function  of  concentration 
by  a  listing  of  ordered  triples  supplied  by  the  user. 

Use  the  following  form! 


first  point 


--  pure  air  y=Q.0>Cc=0.*RH0G=RH0A=  1.0816 


kg/mtt3 


last  point 


-  pure  gas  y=1.0*Cc=RH0E*RH06=RH0E 


A 


■»,M«S75Vv:, 


If  the  release  is  isothermal,  respond  "Y" .  A  positive 
response  causes  DEGAD1SIN  to  ask  for  a  list  of 
concentration,  density,  and  mole  fraction  points  for 
the  gas  mixture.  The  default  response  is  negative. 

For  the  Eagle  series,  the  simulations  were  made  as 
though  the  releases  were  isothermal.  Although  the 
releases  were  not  isothermal,  the  difference  in 
temperature  between  the  liquid  pool  and  the  ambient  air 
were  considered  small  enough  to  be  negligible.  A 
sample  DEGADIS  input  session  for  LNG  can  be  found  in 
Reference  3 . 

Enter  the  three-letter  designation  of  the  diffusing 
gas.  The  properties  of  LNG  as  methane,  LPG  as  propane, 
pure  N02  (without  reaction) ,  and  pure  N204  (without 
reaction)  are  included. 

A  list  of  the  properties  for  the  specified  gas  (if 
available)  is  given.  If  any  of  the  parameters  are  to 
be  changed,  the  first  letter  of  the  parameter  to  be 
changed  in  the  list  is  given  to  the  prompt.  Here,  the 
vapor  density  at  the  storage  temperature  is  changed. 
Note  that  for  N02  (and  N204) ,  the  storage  temperature 
is  assumed  to  be  the  ambient  temperature,  and  the  vapor 
density  is  calculated  by  the  ideal  gas  equation  of 
state.  For  this  simulation,  it  was  unnecessary  to 
change  the  storage  temperature  since  heat  transfer  was 
not  included  although  the  new  vapor  density  value  is 
calculated  from  the  liquid  pool  temperature. 

The  gas  property  list  is  displayed  again.  The  default 
response  is  no  change . 


Since  this  is  an  isothermal  simulation,  the  density  and 
concentration  must  be  specified  as  a  function  of  mole 
fraction. 


Do  you  how  an  input  file  for  the  Density  function?  [a  or  N3 

ENTER  THE  NUMBER  OF  DATA  TRIFLES  (#3^20)  FOR  THE  DENSITY  FUNCTION!  2 

Enter  Hole  frao  Cc  (ha/»tt3) »  then  RHCS  ;kd/mtt3)  ba  triples 
C.»v. <1.0316 

1.»1.?4M.76? 

The  suddested  LOWEST  CONCENTRATION  OF  INTEREST  (aas.lfl/2.) 

:•.?  1.12250E-04  kS/alH.  Enter  the  desired  value!  3.0EH 


Specification  of  source  rate  and  extent. 


Is  this  3  Steady  stste  simulation?  <a  or  N>  Y 


Enter  the  desired  evolution  rate  C=3  kd/sec  !  1.72 
Enter  the  desired  source  radius  [=]  m  !  10. 


In  addition  to  the  inforaation  Just  obtained)  DEGADIS 
reauires  a  series  of  nuaericsl  parameter  files  uhich  use 
the  sane  naae  as  [D1R1RUNNAME  diven  above. 


For  convenience)  example  parameter  files  are  included  for 
each  step .  They  are! 

EXAMPLE. ER1  and 
EXAMPLE. ER2 

Note  that  each  of  these  files  can  be  edited  durind  the  course  of  the 
simulation  if  a  parameter  proves  to  be  out  of  specification. 


Do  you  want 
The  command 


a  command  file  to  be  dene rated  to  execute  the  procedure?  <Y  or  n> 
file  will  be  dene rated  under  the  file  name! 


EAELE6.com 


Do  you  wish  to  initiate  this  procedure?  <y  or  N> 
t 


A  file  may  be  used  to  enter  the  density  information. 
In  the  file,  the  first  line  should  have  the  number  of 
data  triples .  Each  subsequent  line  should  contain  a 
mole  fraction,  concentration,  and  density  for  the  gas 
mixture . 


The  lowest  concentration  of  interest  is  the 
concentration  at  which  the  calculations  are  stopped. 


If  a  steady  state  release  is  to  be  simulated,  type  "Y" 
to  the  prompt.  For  a  steady  simulation,  the  steady 
state  mass  evolution  rate  (kg/s)  and  primary  source 
extent  (m)  are  required. 

A  note  about  the  numerical  parameter  files  is  included. 
These  files  contain  various  constant  values  used  in  the 
programs  to  which  the  user  has  access  without 
recompiling  the  programs.  Access  is  granted  as  a 
convenience . 


DEGADISIN  will  generate  a  command  procedure  suitable 
for  running  the  model  under  VMS. 

If  so  desired,  DEGADISIN  will  initiate  the  command 
procedure  under  VMS.  If  not,  the  program  returns  to 
the  operating  system. 


The  generated  INP  file  for  EAGLE6  is  shown  in  Figure  A-6. 

If  necessary,  the  user  may  edit  the  INP  file  before  beginning 
the  simulation. 

4.  Example  Simulation  Output 

After  proper  completion  of  the  model,  EAGLE6.LIS  contains 
the  output  listing  for  the  steady  state  releases.  Following  is 
a  point-by-point  discussion  of  the  output. 
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UOh_DEGADIS  model  output 


VERSION  1.3 


tmmmmn 

DoLu  input  on 
Source  rrosraa  run  on 


nmttnmm  9-kam986  14:54:06.79  m 

9-MAY-1986  14511:56.90 
9-MAY-1986  14554:06.79 


i 

mmmmmmmmmmmmmmtmmmtmmmmtmm 

t 

%  note: 

*  — 

t 

t  >  All  Calculations  are  liaited  to  circular  lieuid  sources. 

t 

nmmmmmnnmmmtmmmtmtmmmmmmmtnm 


TITLE  block 


EAGLE  6  test  simulation  as  an  ISOTHERMAL  release  of  PURE  H02  WITHOUT 
CHEMICAL  REACTION  usiriS  the  HIM  SOURCE  RATE.  The  source  sas  is  assumed  to 
have  the  sane  tereereture  as  the  LIQUID  POOL  TEMFERATURE  (150. 


Uird  velocity  at  reference  heisht 
Reference  heidiit 


5.59  »/s 

12.00  • 


Surface  roughness  leieth 


1.000E-06  • 


Paiouill  Stability  class 


D 


Mt-nin-Obukhov  lensU. 

infinite 

Cm  i..  .i.-.ii  Jj •.!.!■  ibv!  ■  <’i.  ;on’_t juts- 

Della 

0.13000 

Beta 

0.90000 

Ur  <:  velocity  power  lew  constant 

Alpha 

0.06869 

Tri : l ; o. •  -•Ijcitu 

0.11981 

Ar-bierri  Teneev-*tu.v 
Ambient  Pressure 
Ambient  Absolute  Huaidits 
Vnbient  Relative  Hunidits 


295,75  K 
0.909  ate 
6.758E-03  kS/kS  BDA 
35.00  Z 
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The  date  and  time  DEGADISIN  was  run  is  included. 


The  input  information  gathered  by  DEGADISIN  is  repeated 
to  assist  in  documentation  of  the  simulations. 

Included  here  are  the  Title  Block  and  the  atmospheric 
conditions . 


input : 


© 


Mole  fraction 

0.00000 

1.00000 


CONCENTRATION  OF  t 

y,s/»m 

0.00000 

1,76900 


6AS  DENSITY 
ki/*U2 
1,08160 
1.76900 


Specified  Gas  Properties: 


Molecular  weight! 

46.000 

Storage  temperature: 

295.75  K 

Density  at  storage  temperature  and  ambient  pressure! 

1.7690  kg/att3 

hern  heat  capacity  constant: 

3845.0 

Mean  heat  capacity  power! 

1.0000 

Upper  aole  fraction  contour! 

1.00000E-03 

Lower  *ole  fraction  contour! 

5.00000E-04 

Height  for  isopleths! 

0.50000  a 

Source  input  data  points 


Initial  aass  in  cloud:  0.00000E+00 


TINE 

s 

O.OOOOOE+OO 

6023.0 

6026.0 

6025,0 


SOURCE  STRENGTH 
kS/S 
1.7200 
1.7200 
O.OOOOOE+OO 
O.OOOOOE+OO 


SOURCE  RADIUS 

a 

10.000 

10.000 

O.OOOOOE+OO 

O.OOOOOE+OO 


Calculation  procedure  for  ALPHA*  1 


Entrainaent  prescription  for  PHI!  3 


Layer  thickness  ratio  used  for  average  depth!  2.1500 


Air  entrainment  coefficient  usee  I  0.590 
Gravity  sluitrins  velocity  coefficient  used:  1.150 
Tvjtlieraal  raiculaticn 


Meet  transfer  not  included 


I'.'.IM-  ti ft.*  ini  include1’ 
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Continuing  with  the  input  information,  the  contaminant 
gas  properties  are  output. 


The  specification  of  the  mass  evolution  rate  and  source 
radius  are  output.  For  a  steady- state  release,  there 
is  no  initial  mass  in  the  cloud,  and  the  source 
strength  and  source  radius  are  held  constant  for  an 
arbitrarily  large  period  of  time. 


Finally,  certain  numerical  parameters  and  calculation 
flags  are  displayed.  Some  of  these  are  set  in 
DEGADIS1N  while  others  are  set  in  the  numerical 
parameter  files. 


^T7V.rr/"t 


to  y 

a) 

y  y 

u 

y  x: 

cd 

U  4-1 

y 

y  o 

X 

y 

CO 

o 

E 

3 

y  x> 

c 

rH 

x:  a  • 

O  O 

>> 

y  C  y 

<A4 

G  y  to  3  r-N 

•H  C 

u 

3  3  -h 

o  ta  c  • 

AJ 

cd 

O.  a 

CO 

•rt  -  bO  0 

ca  a) 

co 

co 

y  ^  C  ta 

r-A  A4 

1 

y  3  M 

cd 

3  y  o 

3  Ctf 

3  <d 

y  -h  c 

E  X3  •  x:  to  xi 

O 

w 

>4  3  0 

■rt  N  y  3  y 

r-A  CO 

rj 

3  3  r-t 

<D 

y  \  y  .r-t 

03  Pi 

<  i— 1 

O  y  3 

O 

y  J  -  «a  os 

o  a) 

o 

to 

Pi 

tfl  r-s  Xi 

AJ 

c 

to  y 

2 

•rj  1!  o  y 

a>  a) 

:*.  «  s 

o 

T3  -h  >, 

o  s 

<U 

y  bOi-t 

Ui 

_  X  y  oj  y 

Pi  cd 

B  w 

<t)  y 

c  ^  ta  -h 

3  Pi 

•r4 

•3  y 

<s> 

O  M  y  y  y 

O  cd 

AJ 

G  O  rC 

■h  to  y  x  o 

co  a 

<D 

a  y  o 

u 

y  >-<  y  r-t 

^  O 

y  3 

ta  y  y 

Q)  T> 

O  Pi 

y  o  to 

bO 

y  »h3  > 

rC  Q) 

3 

to  to  c 

£ 

y  w>  o  c 

AJ  AJ 

C  o 

_  -  O 

%r4 

C3S  «  M 

cd 

O  CO 

■O  ^<r-> 

*2 

y  y  -w-  c 

~  rH 

•H 

y  y  y 

£ 

y  -  .y 

a>  3 

P>  >*> 

y  3  y 

rH 

C  T3  C  ^  T3 

CO  o 

U  Pi 

ta  n  c 

O 

o  G  o  >s  ta 

cd  ,-i 

C  (0 

r-t  C  3 

£ 

y  .y  >h  y  y 

a)  cd 

3  TJ 

3  O  ty 

•H 

a  y  .h  y 

r-l  u  ty  C 

y  y 

r-l  C  y  to  o. 

Q) 

o 

r-t  y  0) 

co 

ta  a  <a  c  to 

U  0) 

cd  o 

ca  to  ta 

j-i 

y  o  y  y 

pC 

a> 

y 

0) 

■h  3  y  q  3 

a>  aj 

co  co 

y  y 

AJ 

y  >_r  3 

AJ 

id 

y  x:  3 

0) 

y  y  y  o 

Cd  V4 

0) 

.OHO. 

g 

#XH  ^r-l 

AJ  Q) 

bO  aj 

y  y 

<3 

>  y  o  y  u 

(0  AJ 

£  cd 

,  3 

PA 

S  "rt 

•  Cf4 

•H  AJ 

<y  •  o 

<« 

y  y  to  y 

>*>  cd 

bD  10 

O  T3 

O..C  <8  y  C  X! 

T3 

C  • 

y  y 

y  c  y  y 

Cd 

td  >s 

>•>•3  y 

y  ta  •o 

y  y  x:  T3 

u  3  ta 

V4 

-  y  c  c 

y  y 

o  cd 

ta  r-t 

Ctf 

r-v  y  t-i  y  o 

3  3 

a> 

s  y  y 

•o 

y  y  a  y 

c 

Pi  AJ 

9  ^ 

c 

3  i  3  3  33 

3  'rt 

0)  U) 

3  vt  to 

o 

y  ta  y  y  y 

5 

b0 

to  -H 

u 

w  y  c  x  to 

y  y 

£  0) 

to  y 

<D 

O'  3  o  >y  3 

o  y 

o  ,£ 

<  i-t  X! 

CO 

— '  o.  o  a  xa 

ft.  u 

r-l  AJ 

©  © 


80 


«:.  i 


*0  *  m 

H  O  • 


O  O  II 
•*->  a  n 

CM  — 

jc  o  a 

♦>  I  ro  ^ 

T5  U  J 

^  o 


m  c-i 
«r  m 


«r  00 

*0  “O 


W  N  O' 

ro  -0  O' 


im  -p  m 

O'  o 
•H  M  CM 

O  o  o 


co  >0  -«r 
®  CD  CD 


*<»  Jp  'O 

O'  J>  o- 
ci  cm  cm 


03  H  N 
IO  -0  >0 


m  o  -o 
rs  00  cd 


O'  r--v  -o 

IO  O'  'O 


ki  in  r1 

ro  c- 
c-4  cm  r-j 


O  rH  O 
O  CO  C-4 


O  -O  "O 
O-  CK  O' 
C-4  CM  C-l 


cm  o  im 
rC  co  00 


O'  o  ^ 
V~*  r-4  CM 


*  •  o 

CO  CK  r-l 


r>.  o  cm 
in  o 

N  ^  «f 


O'  o  o 
fH  f*4  CM 


OOO 
C-4  CM  C  4 


co  C' 

*H  ^ 

C-4  CM  C-4 


w'j  tn  p 
CM  pi  CN 


t*i  co  cm 

pm  rn  10 


in  im  O' 
000 
in  -o  n 

000 


~o  ;p  >p 
O'  O'  O' 

CM  c-4  cm 


-O  -O  "O 
O'  O  O 
CM  C  v|  C-4 


05  O  rM 
in  'O  -© 


O'  CO  C' 
O  CM  ^ 
00  co  00 


o  o  o> 

000 

M-  -V  -H 
ui  UJ  Ul 
000 
o  O  o 
o  o  o 


-o  'O  *o 
O'  O'  O' 
CM  C-4  Cl 


m  *-«  cm 
00  'p 

CO  O'  O' 


O  o  0 
O  O  $ 


-O  NO  'p 
O'  t>-  o*. 
C-4  £4  C-4 


111.  6.869E-0J  1.185E-02  1.09  296.  0.000E+00  1.01  18.2  26.7  21.9 
i24.  S.787E-03  9.985E-03  1.09  296.  0.000E+00  1.12  19.6  27.7  22.4 
137.  4.929F-03  9.522E-03  1.08  296.  O.OOOEaOO  1,22  20.9  28,6  22.7 


AD-A172  956 


UNCLASSIFIED 


DEVELOPHENT  OF  VAPOR  DISPERSION  MODELS  FOR  NONNEUTRALLV  2/ 
BUOYANT  GAS  NIXTU.  .  <U>  ARKANSAS  UNIV  FAYETTEVILLE  DEPT 
OF  CHEMICAL  ENGINEERING  T  0  SPICER  ET  AL.  SEP  86 
AFESC/ESL-TR-86-24  DTCG22-88-C-2829  F7G  471  NL 


’St 

> 

'C 

| 

I 


8 


ai 

to 

A 

X 

>v 

0) 

*• 

4) 

41 

r-l 

X 

<A 

•  4) 

0) 

4-1 

•d 

to  <o 

u 

<d 

c 

fsc 

41 

w 

0) 

X  r-l 

•H 

•l-t 

to  41 

td 

d 

> 

•d 

41  O 

W 

td 

4J 

0 

3 

s 

c 

<U 

4> 

■o  g 

to 

cd 

a> 

1-4 

r-4 

e 

<0  A 

C 

c  o 

i— t 

S  X 

0) 

o 

•1-1  <D 

•o 

4J 

o 

<d  4t 

<0 

a> 

41  X 

c 

A 

o 

44 

> 

H 

•w 

C  41 

4> 

4-j 

- 

<0 

3-i-t 

to 

<p 

44 

o 

00 

NX 

1— < 

•H 

a 

10 

>  c 

3 

tO  41 

to 

AJ 

0 

•H 

•rl  *H 

41 

<u 

X 

•o 

p 

•it 

c 

44 

00  « 

X 

c 

60 

-  *f-t 

o 

3 

■d 

C 

41 

•P^ 

•H 

3 

^"v  !* 

•H 

4-1 

o 

w  C 

0 

o 

X 

4-1 

3 

to 

c 

•H 

•rt  o 

It 

44 

41  0) 

td 

O. 

to 

o 

4J 

o 

- 

<0 

X 

■d  x 

4t 

4J 

3 

td 

<0 

c 

41 

41 

•i-t  41 

41 

3 

C 

41 

f— 4 

U  U"l 

o 

c 

& 

e 

o 

i-t 

<d 

3 

u 

•i-t 

® 

VO 

c 

to 

4-1 

44 

* 

O 

3  42 

4-1 

u 

<41  <H 

u 

ai 

e 

41 

4) 

r-H 

O  60 

« 

to 

i— t  <d 

c 

X 

o 

P 

to 

td 

(0  3 

44 

to 

c 

td  4> 

o 

44 

o 

to 

0) 

o 

O 

4-1 

A 

x  P 

a 

o 

44 

0)  44 

e 

4) 

3 

^  o 

41 

e 

(0 

<0 

A  A 

<0 

i — t 

o 

•d 

td 

r-l 

o 

41 

A 

41  41 

o 

C 

o 

X 

to 

X 

r-l 

o 

e 

u 

C 

o 

u 

o 

41 

4J 

<d 

i4 

4-1  <4 

o 

to  i—4 

td 

3  1-4 

4-1 

O 

u 

o 

>4 

o 

« 

44 

a) 

<44 

o 

to 

4l 

• 

(0  •O 

i-t 

41 

U 

> 

O 

•d  C 

41 

3 

f — t 

4)  c 

•d 

O 

id 

0) 

C 

c  s 

c 

4> 

to 

<0  <d 

p 

X 

r-l 

P 

o 

•H  3 

Id 

td 

> 

i 

-H 

41 

i 

o 

i-t 

U  r“l 

c 

4l 

<0 

3  O' 

3 

d  i4 

4J 

£  o 

•1-t 

10 

1-4 

4i 

<0 

• 

(X  P 

4-> 

u 

5  o 

§ 

a 

td  to 

X 

4* 

44 

3 

cd 

o 

0 

3 

s 

•a 

a  c 

41 

X 

3 

o 

V4 

O.T5 

u 

0) 

c 

60 

O 

4t 

V 

• 

c 

41 

3 

(0  3 

o 

1-1 

bO  C 

•o 

®  C 

o 

o 

O.  r-t 

41 

« 

0) 

c 

0  g 

V 

<0 

44 

<d  o 

X 

■ 

•> 

o 

•H 

C  3 

u 

W>rC  u 

«> 

■d 

a> 

c 

g 

<d  r-i 

- 

3 

w 

P 

•d 

« 

c 

0 

P 

41  o 

c 

41 

4> 

i4 

« 

4> 

•H 

o 

to  o 

o 

X 

td 

4t  Xr-4 

41 

<d  i-t 

o 

•rl 

•f-l 

*H 

3  --I 

44 

td 

C 

H 

o 

•d 

•d  u 

41 

e 

■d 

O  r-t 

<0 

y 

•rl 

0) 

(0 

to 

O 

3 

4J  td 

4-1 

•H 

0 

4J 

a> 

V 

4)  4t 

td 

•o 

o 

C  C 

P 

■d 

U 

c 

x 

A  -H 

44 

C  r4 

O  t4 

V 

c 

0) 

<D 

o 

H 

H  4h 

4-t 

o 

U  (4 

0 

•H 

u 

O 

rH 

4.  Model  Limitations  and  Cautions 

DEGADIS  model  application  should  be  limited  to  the 
description  of  atmospheric  dispersion  of  heavier -than- air  gas 
releases  at  ground  level  onto  flat,  unobstructed  terrain  or 
water.  Application  to  releases  from  sources  above  ground  level 
(e.g.  overflow  from  dikes)  would  be  expected  to  give 
conservative  predictions  of  the  downwind  hazard  zones,  but  this 
has  not  been  verified. 

The  dispersion  of  a  heavier- than- air  gas  by  the  action  of 
the  wind  assumes  the  maintenance  of  a  wind  velocity  profile  in 
the  gas  cloud  or  plume  whose  characteristics  are  determined  by 
the  approach  wind  flow  (upwind  of  the  release) .  The  treatment 
of  vertical  momentum  transfer  invokes  the  assumption  of  a 
logarithmic  vertical  velocity  profile,  which  is  in  turn  curve - 
fitted  to  a  power  law  vertical  velocity  profile.  DEGADIS  also 
uses  similarity  forms  for  the  vertical  profile  of  gas 
concentration  in  the  cloud,  and  the  vertical  profile  is 
dependent  on  the  power  law  exponent  o  used  in  the  representation 
of  the  velocity  profile.  The  vertical  velocity  profile,  which 
is  directly  related  to  the  air  entrainment  velocity  into  the 
cloud,  is  dependent  on  the  factors  which  determine  the  structure 
of  the  atmospheric  boundary  surface  layer,  wind  speed,  surface 
roughness,  and  atmospheric  stability.  Consequently,  the 
representations  of  the  vertical  velocity  and  concentration 
profiles  in  DEGADIS  are  subject  to  similar  limitations  as  in 
other  descriptions  of  the  surface  layer.  Table  A- 2  indicates 
typical  recommended  surface  roughness  values.  Table  A-l 
indicates  logarithmic  wind  velocity  profile  corrections  for 
different  atmospheric  stabilities,  along  with  typical  values  of 
the  wind  profile  power  law  exponent  a  determined  in  DEGADIS. 
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